Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

Ok, I will do the same. On Thursday, February 29, 2024 at 12:11:56 AM UTC+5:30 Dima Pasechnik wrote: > On Wed, Feb 28, 2024 at 5:42 PM 'Animesh Shree' via sage-devel < > sage-...@googlegroups.com> wrote: > >> reason scipy factors only square sparse matrices >> >> Problem is basically in

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

On Wed, Feb 28, 2024 at 5:42 PM 'Animesh Shree' via sage-devel < sage-devel@googlegroups.com> wrote: > reason scipy factors only square sparse matrices > > Problem is basically in _superlu.gstrf >

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

reason scipy factors only square sparse matrices Problem is basically in _superlu.gstrf it accepts only one dimension as input rather than both row and col. In c

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

I am currently working for RDF and CDF mainly. I was planning to overload LU function of *sage/matrix/matrix2.pyx* by defining new LU in *sage/matrix/matrix_double_sparce.pyx* . I couldn't understand what is multiprecision version but I guess you are asking for RR or other inexact field.

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

As Dima points out, SPLU decomposition is probably more useful in practice, so it's great to expose it. It should probably be exposed through a different method than PLU decomposition, though. For the return format: permutations are much more compact and efficient to describe permutations than

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

Yes Because of same reason I tried to commented scipy code to test this. I got some error saying *RuntimeError: Factor is exactly singular* But same worked for

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

There is a good reason for numerics people to adopt "SuperLU" factorisations over the classical PLU sparse decomposition - namely, it's more stable. Perhaps it should be made the Sage's default for sparse RDF matrices, too. By the way, https://portal.nersc.gov/project/sparse/superlu/superlu_ug.pdf

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

I went through the link. It also returns perm_c and perm_r and the solution is represented as Pr * (R^-1) * A * Pc = L * U It is similar to one returned by scipy >>> lu.perm_r array([0, 2, 1, 3], dtype=int32) >>> lu.perm_c array([2, 0, 1, 3], dtype=int32) I think it doesn't

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

One more option would be umfack via scikits.umfpack: https://scikit-umfpack.github.io/scikit-umfpack/reference/scikits.umfpack.UmfpackLU.html Regards, Max On Wednesday, February 28, 2024 at 7:07:53 AM UTC-5 Animesh Shree wrote: > One thing I would like to suggest. > > We can provide multiple

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

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 (Note - link

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

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

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

On 27 February 2024 17:25:51 GMT, 'Animesh Shree' via sage-devel 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

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

See the documentation: https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.SuperLU.html As you can see there it computes permutation matrices Pr and Pc and lower and upper triangular matrices L,U such that Pr*A*Pc = L*U so it's not computing a normal LU decomposition of a

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

Sorry for multiple messages I just want to say >sage: p,l,u = a.LU(force=True) >sage: p >{'perm_r': [1, 0], 'perm_c': [1, 0]} It ( {'perm_r': [1, 0], 'perm_c': [1, 0]} ) represents transpose and it cannot be represented as permutation matrix. Similar cases may arise for other matrices. On

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

For transpose : In the example we can see permutations are provided as arrays for rows and cols. The permutation is equivalent of taking transpose of matrix. But we cant represent transpose as a permutation matrix. >>> a = np.matrix([[1,2],[3,5]]) >>> # a * perm = a.T >>> # perm =

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

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

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

On 27 February 2024 15:34:20 GMT, 'Animesh Shree' via sage-devel 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,

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

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:

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

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 0andU11 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

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

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

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

Le 05/02/2024 à 11:30, Dima Pasechnik a écrit : A good implementation of LU decomposition ought actually to take non-square matrix as input, and have the indices adjusted appropriately in the algorithm, so it's indeed a bit strange that superLU only takes square matrices (?). Perhaps

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

On 5 February 2024 06:30:45 GMT, 'Animesh Shree' via sage-devel wrote: >The other library that scipy uses is SuperLU : >https://portal.nersc.gov/project/sparse/superlu/ >for the function scipy.sparse.linalg.splu :

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

On 5 February 2024 05:36:45 GMT, 'Animesh Shree' via sage-devel wrote: >I tried to use ctypes pythonapi to extract PyCapsule_New objects using >PyCapsule_GetPointer. But got erorr : > *ArgumentError: argument 1: TypeError: Don't know how to >convert parameter 1* >Attaching

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

I tried to use ctypes pythonapi to extract PyCapsule_New objects using PyCapsule_GetPointer. But got erorr : *ArgumentError: argument 1: TypeError: Don't know how to convert parameter 1* Attaching file

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

The other library that scipy uses is SuperLU : https://portal.nersc.gov/project/sparse/superlu/ for the function scipy.sparse.linalg.splu : https://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.splu.html#scipy-sparse-linalg-splu Scipy supports only factorization for square

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

IMHO opaque C objects may be accessed from Python, it just needs an implementation On 29 January 2024 22:51:48 GMT, Nils Bruin wrote: >By the looks of it, the routines you'd be using would be coming from >umfpack. cvxopt has chosen to package the details of LU factorization in >opaque

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

By the looks of it, the routines you'd be using would be coming from umfpack. cvxopt has chosen to package the details of LU factorization in opaque objects and instead offer routines to use these decompositions (via taking the opaque object as input). In scipy, umfpack is also used:

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

It returns PyCapsule object, I am looking if we can decapsulate it in python/cython. If it doesn't work then I will go for Scipy On Friday, January 19, 2024 at 4:04:37 PM UTC+5:30 Animesh Shree wrote: > Thank You I am looking into it > > On Friday, January 19, 2024 at 3:43:29 PM UTC+5:30 Dima

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

I tried using that but I am unable to access the factors. I asked them( CVXOPT ) about this. For now I am falling back to Scipy implementation. But it ( scipy.sparse.linalg.splu

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

Thank You I am looking into it On Friday, January 19, 2024 at 3:43:29 PM UTC+5:30 Dima Pasechnik wrote: > On Fri, Jan 19, 2024 at 9:35 AM 'Animesh Shree' via sage-devel > wrote: > > > > I have been looking into sparse matrix LU decomposition issue. > > > > In that conversation SuiteSparse is

Re: [sage-devel] SuiteSparse and sage and sparse_matrix.LU()

On Fri, Jan 19, 2024 at 9:35 AM 'Animesh Shree' via sage-devel wrote: > > I have been looking into sparse matrix LU decomposition issue. > > In that conversation SuiteSparse is proposed to handle sparse matrix > decomposition. > > What I could find that it is mostly C and Matlab. I don't know

[sage-devel] SuiteSparse and sage and sparse_matrix.LU()

I have been looking into sparse matrix LU decomposition issue . In that conversation SuiteSparse is proposed to handle sparse matrix decomposition. What I could find that it is