From e1fc87a96411d2364e1b5f1bda632b297a9ea537 Mon Sep 17 00:00:00 2001 From: martinjrobins Date: Thu, 3 Mar 2022 17:14:44 +0000 Subject: [PATCH] #1863 fix scikits solvers --- pybamm/solvers/algebraic_solver.py | 3 ++- pybamm/solvers/scikits_ode_solver.py | 8 ++++++-- tests/unit/test_solvers/test_scikits_solvers.py | 12 ++++++------ 3 files changed, 14 insertions(+), 9 deletions(-) diff --git a/pybamm/solvers/algebraic_solver.py b/pybamm/solvers/algebraic_solver.py index 1aa08640a7..81f0348c6f 100644 --- a/pybamm/solvers/algebraic_solver.py +++ b/pybamm/solvers/algebraic_solver.py @@ -67,7 +67,8 @@ def _integrate(self, model, t_eval, inputs_dict=None): y0 = model.y0 if isinstance(y0, casadi.DM): - y0 = y0.full().flatten() + y0 = y0.full() + y0 = y0.flatten() # The casadi algebraic solver can read rhs equations, but leaves them unchanged # i.e. the part of the solution vector that corresponds to the differential diff --git a/pybamm/solvers/scikits_ode_solver.py b/pybamm/solvers/scikits_ode_solver.py index 847d44badc..addeb728f5 100644 --- a/pybamm/solvers/scikits_ode_solver.py +++ b/pybamm/solvers/scikits_ode_solver.py @@ -98,8 +98,12 @@ def _integrate(self, model, t_eval, inputs_dict=None): events = model.terminate_events_eval jacobian = model.jac_rhs_eval - def eqsydot(t, y, return_ydot): - return_ydot[:] = derivs(t, y, inputs) + if model.convert_to_format == "casadi": + def eqsydot(t, y, return_ydot): + return_ydot[:] = derivs(t, y, inputs).full().flatten() + else: + def eqsydot(t, y, return_ydot): + return_ydot[:] = derivs(t, y, inputs).flatten() def rootfn(t, y, return_root): return_root[:] = [event(t, y, inputs) for event in events] diff --git a/tests/unit/test_solvers/test_scikits_solvers.py b/tests/unit/test_solvers/test_scikits_solvers.py index de4f496ccb..c4bda26a7f 100644 --- a/tests/unit/test_solvers/test_scikits_solvers.py +++ b/tests/unit/test_solvers/test_scikits_solvers.py @@ -49,10 +49,10 @@ class Model: length_scales = {} convert_to_format = "python" - def residuals_eval(self, t, y, ydot, inputs): - return np.array([0.5 * np.ones_like(y[0]) - ydot[0], 2 * y[0] - y[1]]) + def rhs_algebraic_eval(self, t, y, inputs): + return np.array([0.5 * np.ones_like(y[0]), 2 * y[0] - y[1]]) - def jacobian_eval(self, t, y, inputs): + def jac_rhs_algebraic_eval(self, t, y, inputs): return np.array([[0.0, 0.0], [2.0, -1.0]]) model = Model() @@ -101,12 +101,12 @@ class Model: convert_to_format = "python" len_rhs_and_alg = 2 - def residuals_eval(self, t, y, ydot, inputs): + def rhs_algebraic_eval(self, t, y, inputs): return np.array( - [0.5 * np.ones_like(y[0]) - 4 * ydot[0], 2.0 * y[0] - y[1]] + [0.5 * np.ones_like(y[0]), 2.0 * y[0] - y[1]] ) - def jacobian_eval(self, t, y, inputs): + def jac_rhs_algebraic_eval(self, t, y, inputs): return np.array([[0.0, 0.0], [2.0, -1.0]]) model = Model()