| Topic |
Links |
|---|---|
| Code and Changes | Submission via Patch Tracker |
| Community and Discussion | Discourse Topic - Ode15{i,s} Project Page |
| Repository | Octave Github Mirror |
My project was Ode15{i,s} : Matlab Compatible Solvers, and as the name suggests revolved around two specific solvers from the set provided by Octave for initial value problems for Ordinary Differential Equations that have a MATLAB-compatible interface. My project was a continuation to the initial implementation of these during GSoC ‘16, and focused on improving the current state of the solvers, heavily involved with the SUNDIALS Library used in the backend.
The main problem to be solved stemmed from the fact that SUNDIALS IDA made use of the KLU Library for sparse computations. Any issue in linking this with respect to SUNDIALS, so as to generate the SUNLINSOL_KLU Linear Solver Interface, which packaged versions of SUNDIALS are often not compiled with (yes, even for widely popular operating systems such as Ubuntu 20.04), this resulted in the entire computation being carried out using dense Jacobian matrices instead, causing the user huge hits in terms of performance !
Hence, this project effectively created an interface to SUNDIALS IDA API, using Octave’s own classes and solvers, also simultaeneously taking advantage of not having to create copies of the data in SUNDIALS as well as Octave’s vectorized operations, creating a safety net for the end user such that sparse computations are never disabled in the Ode15{i,s} solvers.
As is for most users that use a user interface, what happens in the backend, stays in the backend.
Meaning, before the existence of my project, if the KLU library was not found, the user would merely recieve the warning
KLU library not found. This will result in some lack of functionality for sparse matrices.
What this lack of functionality actually meant was, internally, even when the user entered a “sparse” Jacobian, it was stored and computed using a “dense representation”, causing extreme performance degradation.
Now, the warning has been replaced with
KLU library not found. This will result in some lack of functionality for sparse matrices.
Sparse computations in ode15{i,s} solvers will still be supported but will suffer a hit in time taken.
Here, ode15{i,s} still supports “sparse computations”, meaning the user-entered Jacobian is internally stored in a “CSC (Compressed Sparse Column)” representation, and solved using Octave’s own polymorphic sparse LU solver, instead of KLU.
So, to what extent does this help in the performance degradation taking place earlier?
Solving a PDE with variable NDOFS (3D Heat Equation) using Ode15i
The test file used can be found
here
KLU Solversindicates SUNDIALS IDA using the SUNLinsol_KLU module as the linear solver, which is an interface to Suitesparse/KLU
__________________________________
Octave Solversindicates SUNDIALS IDA using my custom-made OCTLinSol_Gen module as the linear solver, an interface to Octave’s polymorphic solver specified above, which is now used for sparse computation if KLU is not found.
__________________________________
Dense Solversindicates SUNDIALS IDA using the SUNLinSol_Dense module as the linear solver, which was the default earlier if KLU was not found.
N = 500 ndof = 343
KLU Solvers : Elapsed time is 0.0539529 seconds.
Octave's Solvers : Elapsed time is 0.124162 seconds.
Dense Solvers : Elapsed time is 0.389471 seconds.
N = 1000 ndof = 729
KLU Solvers : Elapsed time is 0.183311 seconds.
Octave's Solvers : Elapsed time is 1.42844 seconds.
Dense Solvers : Elapsed time is 3.63896 seconds.
N = 5000 ndof = 4913
KLU Solvers : Elapsed time is 8.73904 seconds.
Octave's Solvers : Elapsed time is 120.542 seconds.
Dense Solvers : Elapsed time is 449.841 seconds.
N = 10,000 ndof = 9261
KLU Solvers : Elapsed time is 62.2595 seconds.
Octave's Solvers : Elapsed time is 182.939 seconds.
Dense Solvers : Elapsed time is 2189.79 seconds.
These results show that while KLU is still better, and is hence the first choice, our implementation of using the Sparse LU solver is hugely beneficial when KLU is not found, preventing the solver from falling to the dense implementation.
Learnt about creating custom implementations of required data structures to interface with SUNDIALS IDA API.
Created custom implementations for - NVector, SUNMatrix and SUNLinearSolver as wrappers to corresponding classes - ColumnVector, SparseMatrix and SparseMatrix->solve() in Octave. Added function-wise tests to test each implementation individually, and successfully tested interface with IDA API via the IDA2DHeatKry example.
The custom implementations created have now been consolidated and stored in Octave’s file structure as
libinterp/dldfcn/oct-sundials.h
libinterp/dldfcn/oct-sundials.cc
Integrated custom implementations and tested the performance of interface to Octave’s Solver against KLU in Ode15i. Integrated with current build system that uses GNU Autotools, compiling __ode15__.cc oct-file linked against custom implementations.
On obtaining the shown results and deciding to use Octave’s Solvers as default only when KLU is not found, made use of particular header check macros set by autoheader, to redesign the oct-file as follows :
#if defined (HAVE_SUNDIALS)
#include "oct-sundials.h"
/* containing NVector_Octave, OCTMatrix_Sparse
and OCTLinsol_Sparse */
if (user passes sparse jacobian)
{
#if defined (HAVE_SUNDIALS_SUNLINSOL_KLU)
use
- NVector_Serial
- SUNMatrix_Sparse
and
- SUNLinSol_KLU
#else
use
- NVector_Octave
- OCTMatrix_Sparse
and
- OCTLinSol_Sparse
#endif
}
else (user has passed dense jacobian)
{
use
- NVector_Dense
- SUNMatrix_Dense
and
- SUNLinSol_Dense
}
The custom implementations and their interface to Ode15{i,s} solvers is hence complete, such that any user could use this by merely running configure and make, and it also passes all the previous tests in the ode15{i,s} m-files with the new implementations, keeping KLU optional and using it when required.
Some improvements that could be made are :
For example, the #HAVE_SUNDIALS macro currently checks for a lot of Sundials interfaces, such as NVecSerial, which might not be needed if we are to just use the custom implementations, which only require basic Sundials files. Simpler macros code be made for this usecase.
after which the code should prove to be a useful addition to the Matlab-Compatible ODE and DAE Solvers offered in Octave.
After my exhilirating journey of
GSoC @ GNU Octavein the past few months, I am proud of a lot of things, but mostly of myself, for having been able to absorb even tiny bits from the vast knowledge of my mentors and peers, this priceless experience of each up, down, error and success in the lines of code written by me giving me much to look forward to in a future career in development, research and open-source.
I am grateful and extremely thankful to the entire community at GNU Octave, with everyone from experienced Octave developers to even fellow GSoC contributors never hesitant to contribute to every discussion, ranging from mathematical code at the heart of Octave, to build systems and debugging errors. I have carried out most of the discussion throughout my project at the developer forum and am truly amazed at the inclusivity and collaborative spirit of everyone here, making my first deep-dive into open-source an extremely cherished experience!
But most of all I would like to thank my mentor, Carlo de Falco for his expert guidance and support. My over-ambitiousness and naivety of including too much in my intial proposal were met with calm advice to be accurate and even conservative when estimating the effort of tasks so as to be able to work more efficiently, my questions were met with patient and detailed explanations, not just of the answer at hand but also of how to be able to arrive at it myself the next time, and most importantly, my mistakes were met with concise criticism rather than reprimandation.
Instead of hard-set deadlines to meet, I was encouraged and set achievable goals for myself, even being met with content when I redefined the timeline in the first few weeks of the project. This allowed me to truly explore, learn and develop a sense of accountability to complete my project, being allowed to learn myself about simple things such as name-mangling by C++ compilers or artificial polymorphism in C to the importance of vectorized operations and memory management (accompanied with a good many nights spent with GDB), and as a result I have improved my programming skills significantly since when I started, learning both theory and by practice !
I hope to see this project to its completion and keep contributing to GNU Octave.