IDAKLU Solver#

class pybamm.IDAKLUSolver(rtol=0.0001, atol=1e-06, root_method='casadi', root_tol=1e-06, extrap_tol=None, output_variables=None, options=None)[source]#

Solve a discretised model, using sundials with the KLU sparse linear solver.

Parameters:
  • rtol (float, optional) – The relative tolerance for the solver (default is 1e-6).

  • atol (float, optional) – The absolute tolerance for the solver (default is 1e-4).

  • root_method (str or pybamm algebraic solver class, optional) – The method to use to find initial conditions (for DAE solvers). If a solver class, must be an algebraic solver class. If “casadi”, the solver uses casadi’s Newton rootfinding algorithm to find initial conditions. Otherwise, the solver uses ‘scipy.optimize.root’ with method specified by ‘root_method’ (e.g. “lm”, “hybr”, …)

  • root_tol (float, optional) – The tolerance for the initial-condition solver (default is 1e-6).

  • extrap_tol (float, optional) – The tolerance to assert whether extrapolation occurs or not (default is 0).

  • output_variables (list[str], optional) – List of variables to calculate and return. If none are specified then the complete state vector is returned (can be very large) (default is [])

  • options (dict, optional) –

    Addititional options to pass to the solver, by default:

    options = {
        # Print statistics of the solver after every solve
        "print_stats": False,
        # Number of threads available for OpenMP
        "num_threads": 1,
        # Evaluation engine to use for jax, can be 'jax'(native) or 'iree'
        "jax_evaluator": "jax",
        ## Linear solver interface
        # name of sundials linear solver to use options are: "SUNLinSol_KLU",
        # "SUNLinSol_Dense", "SUNLinSol_Band", "SUNLinSol_SPBCGS",
        # "SUNLinSol_SPFGMR", "SUNLinSol_SPGMR", "SUNLinSol_SPTFQMR",
        "linear_solver": "SUNLinSol_KLU",
        # Jacobian form, can be "none", "dense",
        # "banded", "sparse", "matrix-free"
        "jacobian": "sparse",
        # Preconditioner for iterative solvers, can be "none", "BBDP"
        "preconditioner": "BBDP",
        # For iterative linear solver preconditioner, bandwidth of
        # approximate jacobian
        "precon_half_bandwidth": 5,
        # For iterative linear solver preconditioner, bandwidth of
        # approximate jacobian that is kept
        "precon_half_bandwidth_keep": 5,
        # For iterative linear solvers, max number of iterations
        "linsol_max_iterations": 5,
        # Ratio between linear and nonlinear tolerances
        "epsilon_linear_tolerance": 0.05,
        # Increment factor used in DQ Jacobian-vector product approximation
        "increment_factor": 1.0,
        # Enable or disable linear solution scaling
        "linear_solution_scaling": True,
        ## Main solver
        # Maximum order of the linear multistep method
        "max_order_bdf": 5,
        # Maximum number of steps to be taken by the solver in its attempt to
        # reach the next output time.
        # Note: this value differs from the IDA default of 500
        "max_num_steps": 100000,
        # Initial step size. The solver default is used if this is left at 0.0
        "dt_init": 0.0,
        # Maximum absolute step size. The solver default is used if this is
        # left at 0.0
        "dt_max": 0.0,
        # Maximum number of error test failures in attempting one step
        "max_error_test_failures": 10,
        # Maximum number of nonlinear solver iterations at one step
        # Note: this value differs from the IDA default of 4
        "max_nonlinear_iterations": 40,
        # Maximum number of nonlinear solver convergence failures at one step
        # Note: this value differs from the IDA default of 10
        "max_convergence_failures": 100,
        # Safety factor in the nonlinear convergence test
        "nonlinear_convergence_coefficient": 0.33,
        # Suppress algebraic variables from error test
        "suppress_algebraic_error": False,
        ## Initial conditions calculation
        # Positive constant in the Newton iteration convergence test within the
        # initial condition calculation
        "nonlinear_convergence_coefficient_ic": 0.0033,
        # Maximum number of steps allowed when `init_all_y_ic = False`
        # Note: this value differs from the IDA default of 5
        "max_num_steps_ic": 50,
        # Maximum number of the approximate Jacobian or preconditioner evaluations
        # allowed when the Newton iteration appears to be slowly converging
        # Note: this value differs from the IDA default of 4
        "max_num_jacobians_ic": 40,
        # Maximum number of Newton iterations allowed in any one attempt to solve
        # the initial conditions calculation problem
        # Note: this value differs from the IDA default of 10
        "max_num_iterations_ic": 100,
        # Maximum number of linesearch backtracks allowed in any Newton iteration,
        # when solving the initial conditions calculation problem
        "max_linesearch_backtracks_ic": 100,
        # Turn off linesearch
        "linesearch_off_ic": False,
        # How to calculate the initial conditions.
        # "True": calculate all y0 given ydot0
        # "False": calculate y_alg0 and ydot_diff0 given y_diff0
        "init_all_y_ic": False,
        # Calculate consistent initial conditions
        "calc_ic": True,
    }
    

    Note: These options only have an effect if model.convert_to_format == ‘casadi’

Extends: pybamm.solvers.base_solver.BaseSolver

jaxify(model, t_eval, *, output_variables=None, calculate_sensitivities=True, t_interp=None)[source]#

JAXify the solver object

Creates a JAX expression representing the IDAKLU-wrapped solver object.

Parameters:
  • model (pybamm.BaseModel) – The model to be solved

  • t_eval (numeric type, optional) – The times at which to stop the integration due to a discontinuity in time.

  • output_variables (list of str, optional) – The variables to be returned. If None, all variables in the model are used.

  • calculate_sensitivities (bool, optional) – Whether to calculate sensitivities. Default is True.

  • t_interp (None, list or ndarray, optional) – The times (in seconds) at which to interpolate the solution. Defaults to None, which returns the adaptive time-stepping times.

set_up(model, inputs=None, t_eval=None, ics_only=False)[source]#

Unpack model, perform checks, and calculate jacobian.

Parameters:
  • model (pybamm.BaseModel) – The model whose solution to calculate. Must have attributes rhs and initial_conditions

  • inputs (dict, optional) – Any input parameters to pass to the model when solving

  • t_eval (numeric type, optional) – The times at which to stop the integration due to a discontinuity in time.