Hi Oscar, Thank you for your detailed feedback! I now see that rather than adding multiple matrix decompositions, a well-optimized LU and Fraction-Free LU (FFLU) implementation is far more impactful in sympy. Now that I have gotten a suitable topic that is GSOC worthy and it is worth spending time on during the summer. Here is a draft template of my idea, I await your feedback for refinements.
Title implementing division based LU & Fraction-Free LU (FFLU) algorithms in SymPy Idea LU decomposition is a fundamental operation in linear algebra, used for solving systems of linear equations, computing matrix inverses, and finding determinant values. Currently, SymPy lacks a fully optimized and integrated LU decomposition, especially for fraction-free computations. This project aims to: - Implement both sparse and dense LU/FFLU in DomainMatrix, ensuring efficient computations for various matrix types. - Make use of python-flint for FLINT's FFLU function when possible and ensure that the outputs are comparable to the SymPy implementation. - Implement numerically stable LU decomposition version for floats - Improve Matrix methods that call into the DomainMatrix methods, making LU/FFLU easily accessible to users and improving overall computational speed. - Develop upper/lower triangular solvers for LU/FFLU. - Solve over/under determined systems - Integrate LU-based solving into core SymPy functions (solve, Matrix.solve, linsolve, _linsolve_aug) for better performance. - Utilize LU for computing matrix inverse, RREF, and null space efficiently. - Add Extensive benchmarking & lots of timing to compare performance of different approaches for different domains, shapes, densities etc and will Decide when different algorithms should be used. By completing these improvements, LU/FFLU will become the backbone of SymPy’s equation-solving infrastructure, making linear algebra operations significantly faster and more robust. Status - SymPy currently has basic LU decomposition, but it lacks an optimized fraction-free version and efficient sparse implementations. - Work on QRD using FFLU is in progress (see PR #27423 <https://github.com/sympy/sympy/pull/27423>) - (PR #26000 <https://github.com/sympy/sympy/pull/26000>) is refactoring SymPy’s linear equation solvers. This project will extend and integrate LU/FFLU into that. - Python-Flint provides FFLU for integer matrices, and this project will ensure SymPy uses FLINT’s version when applicable. Involved Software - SymPy - python-flint Difficulty Advanced - Matrix decomposition algorithms (LU, FFLU, QRD). - Symbolic and exact arithmetic (ZZ, QQ, python-flint). - Sparse & dense matrix operations in SymPy’s DomainMatrix. - Integration with existing linear solvers in linsolve, _linsolve_aug, and Matrix.solve. - Optimizing performance for exact arithmetic (ZZ, QQ) vs floating-point numbers (RR, QQ) is non-trivial. Prerequisite Knowledge - knowledge of linear algebra, especially LU decomposition and fraction-free algorithms. - Needs understanding of sparse vs. dense matrix operations. Project Length 175/350hours On Fri, 7 Feb 2025 at 18:14, Oscar Benjamin <oscar.j.benja...@gmail.com> wrote: > Hi Temiloluwa, > > Yes, these kinds of things would be good to have in SymPy. More > important than adding many different factorisations though is having > good implementations of the most important ones. I would say that just > implementing fraction-free and division-based LU algorithms fully > would make a suitable GSOC project. I suspect that you currently > underestimate how much is involved in that even though you have a PR > that implements the algorithm already. The additional things that are > needed to complete the LU implementation are things like: > > - Both sparse and dense implementations of LU and FFLU for > DomainMatrix (consider if any very different sparse algorithms are > worth implementing). > - Make use of python-flint for FLINT's FFLU function when possible and > ensure that the outputs are comparable to the SymPy implementation. > - Numerically stable version of LU for floats. > - Matrix methods that call into the DomainMatrix methods making them > accessible to users and speeding up end user calculations. > - Upper and lower triangular solves for LU and fraction-free LU. > - Solving over/under determined systems. > - Make use of LU to compute other things like matrix inverse, RREF, > nullspace. > - Lots of timings and benchmarks to compare performance of different > approaches for different domains, shapes, densities etc. Decide when > different algorithms should be used. > - Make it so that solve and Matrix.solve and other functions use > LU-based solve if it is better. > > See also this PR which I did not finish but was intended to make all > linear equation solving code go through a single place so that it > could be optimised for all cases: > > https://github.com/sympy/sympy/pull/26000 > > Putting all of that sort of stuff together is more useful than adding > other less frequently used factorisations. A good implementation of LU > is a foundation for many other things so it is worth spending some > time to make it as good as possible. As you have seen other things > like QRD can be almost trivial to implement if they can leverage the > LU decomposition. > > In general not all matrix factorisations that can be defined > abstractly are actually suitable for exact or symbolic arithmetic > rather than floating point arithmetic and vice versa. For example I am > not sure that I have seen any implementation of the polar > decomposition for exact numbers. Unless SymPy could provide something > that is more useful than SciPy's polar function I'm not sure it is > worth adding. Note that a high precision implementation should be part > of mpmath rather than SymPy. The only value implementing this in SymPy > would be if it could be exact but I don't know how you would compute > the polar decomposition exactly or whether that is something that only > makes sense for very small matrices. > > On Fri, 31 Jan 2025 at 21:13, Temiloluwa <ytemilol...@gmail.com> wrote: > > > > Hello SymPy Community, > > > > My name is Temiloluwa, and I would like to propose the addition of more > matrix decomposition methods (Schur, Polar, and Hermite Decomposition) as > my project idea for Google Summer of Code (GSoC) 2025. > > > > I have successfully gotten a PR in for a QR decomposition method for > DomainMatrix and I am currently finalizing a PLDU decomposition and a > fraction-free QR decomposition (QRD) pull request for DomainMatrix, aimed > at improving the GramSchmidt process > > > > From my exploration of the SymPy codebase, I have observed that the few > matrix decomposition methods that exist are LU, QR, Cholesky variants and > some others to name a few. However, important methods like: > > > > Schur Decomposition (A=QTQH) > > Polar Decomposition (A=UP) > > Hermite Decomposition (A=LDLH) > > > > are yet to be implemented. > > > > Integrating these decompositions would enhance SymPy’s linear algebra > capabilities, particularly in eigenvalue computation and numerical > stability. I look forward to discussing this further and receiving feedback > from the community. > > > > Best regards, > > Temiloluwa > > > > -- > > You received this message because you are subscribed to the Google > Groups "sympy" group. > > To unsubscribe from this group and stop receiving emails from it, send > an email to sympy+unsubscr...@googlegroups.com. > > To view this discussion visit > https://groups.google.com/d/msgid/sympy/008ee22b-e3c3-40b0-8565-18e2f1aef17cn%40googlegroups.com > . > > -- > You received this message because you are subscribed to the Google Groups > "sympy" group. > To unsubscribe from this group and stop receiving emails from it, send an > email to sympy+unsubscr...@googlegroups.com. > To view this discussion visit > https://groups.google.com/d/msgid/sympy/CAHVvXxTQX8t9LB-cT5RfOv-t1vDZx%2B9PQD55bTht-gQr33XN8w%40mail.gmail.com > . > -- You received this message because you are subscribed to the Google Groups "sympy" group. To unsubscribe from this group and stop receiving emails from it, send an email to sympy+unsubscr...@googlegroups.com. To view this discussion visit https://groups.google.com/d/msgid/sympy/CAP-ouaXQU%2BD3Z86AhxUV4j7jV49__RkMxxKGjiN2X_%3D8Jbuicw%40mail.gmail.com.
