modelon-community / Assimulo

Assimulo is a simulation package for solving ordinary differential equations.
https://jmodelica.org/assimulo/index.html
GNU Lesser General Public License v3.0
66 stars 16 forks source link

KLU sparse solver and SPGMR preconditioning for IDA solver #86

Open arajaku opened 4 months ago

arajaku commented 4 months ago

I am trying to pass Jacobian as a sparse matrix through a jac function for the IDA solver. However, the code only take two types of linear solvers i.e. dense or SPGMR.

I understand it requires the interface for sparse KLU solver library in IDA which is not yet implemented in Assimulo. So, I proceeded with the SPGMR solver to pass a preconditioned matrix as per SUNDIALS documentation.

In the Assimulo documentation examples, I could find an example (ida_with_jac_spgmr.py) where jacv (Jacobian times vector) function is implemented. In this example, I figured there is no preconditioning used. Then looked at another example in CVODE (Explicit_Problem) where preconditioning is used (cvode_with_preconditioning.py).

However, trying out the same structure for IDA (Implicit_Problem) didn't work. The error was, the options for preconditioning were not available for IDA solver. I came to a conclusion about the behavior of this by looking into the source code of sundials.pyx.

---------------------------------Snippet of source code-------------------------------------------------------------------------

# From [sundials.pyx file](https://github.com/modelon-community/Assimulo/blob/master/src/solvers/sundials.pyx)---provides only two solvers-DENSE/SPGMR
if self.options["linear_solver"] == 'DENSE':

elif self.options["linear_solver"] == 'SPGMR':
    # PREC_NONE implies no preconditioner option is supplied to SUNDIALS
    IF SUNDIALS_VERSION >= (3,0,0):
        #Create the linear solver
        IF SUNDIALS_VERSION >= (4,0,0):
            IF SUNDIALS_VERSION >= (6,0,0):
                self.sun_linearsolver = SUNDIALS.SUNLinSol_SPGMR(self.yTemp, **PREC_NONE**, 0, ctx)
            ELSE:
                self.sun_linearsolver = SUNDIALS.SUNLinSol_SPGMR(self.yTemp, **PREC_NONE**, 0)
        ELSE:
            self.sun_linearsolver = SUNDIALS.SUNSPGMR(self.yTemp, PREC_NONE, 0)

        #Attach it to IDAS
        IF SUNDIALS_VERSION >= (4,0,0):
            flag = SUNDIALS.IDASetLinearSolver(self.ida_mem, self.sun_linearsolver, NULL)
        ELSE:
            flag = SUNDIALS.IDASpilsSetLinearSolver(self.ida_mem, self.sun_linearsolver)

From the comments I added in the code snippet from the source code, I understood that the preconditioning functionality for the SPGMR solver is not exposed in sundials.pyx file.

Let me know if I my understanding is correct. If yes, is there a chance for this functionality to be exposed? It would be really helpful if you could provide some information on this.

PeterMeisrimelModelon commented 4 months ago

Hi,

I agree with your assessment of the current status, IDA + SPGMR does currently not support preconditioning.

There are currently no plans to enhance support for this from our side. There is however the possibility for you to work on this and contribute to the development of Assimulo.

/Peter