[Bug]: IDAKLU SOLVER output_variables long simulation and state issues

[aaltno]

PyBaMM Version

24.11.2

Python Version

3.12.4

Describe the bug

IDAKLU SOLVER output_variables long simulation and state issues:

Hi! I have been running ageing simulations with CASADI for a while and wanted to use IDAKLU for CC ageing or calendar ageing as it is faster. With the new option for output_variables, the memory allocation is much smaller (x10), good for parallel or threaded computing.

Now, two issues arise. Importantly, none of these errors occur with Casadi or IDAKLU (when output_variables are not specified)

So basically using last_state in looped simulations and using .cycles to retrieve cycle specif data does not work too well with IDAKLU when output_variables is active :)

1) When using cycles and selecting output_variables, one cannot use the previous way of working:
for example like this: simulation_IDA3_checkup.solution.cycles[-1]['Throughput capacity [A.h]'].entries
but instead must go with a work around like this:
capacity_full_checkup = simulation_IDA3_checkup.solution.get_data_dict()['Throughput capacity [A.h]']
cycles = simulation_IDA3_checkup.solution.get_data_dict()['Cycle']
capacity_discharge_section= [d for d, c in zip(capacity_full_checkup, cycles) if c == 2]
I find this work around inconvenient and would prefer that what works with CASADI and IDAKLU (without output_Variables selected) works also with output_variables selected

2) when using cycling and for example using last states, I get an error with starting_solution
example:

simulation_IDA3_ageing_cylic.solve(starting_solution=last_state)    
last_state = simulation_IDA3_ageing_cylic.solution.last_state
simulation_IDA3_checkup.solve(starting_solution=last_state)

this will give me an error if simulation_IDA3_ageing_cylic = pybamm.Simulation(... , solver=solver_IDAKLU2 , ....)
but it will not give me an error if simulation_IDA3_ageing_cylic = pybamm.Simulation(... , solver=solver_IDAKLU1 , ....)
solver_IDAKLU1 = pybamm.IDAKLUSolver(rtol=0.00001, atol=1e-07, root_method='casadi', root_tol=1e-07, extrap_tol=None)
solver_IDAKLU2 = pybamm.IDAKLUSolver(rtol=0.00001, atol=1e-07, root_method='casadi', root_tol=1e-07, extrap_tol=None, output_variables=['Time [s]', 'Current [A]', "Voltage [V]" ], options=None)

for this last example I have a work around which is to use starting_solution = simulation_IDA3_ageing_calendar.solution, so using the whole solution instead of the last_state (which is I think more efficient and has a smaller memory footprint)

for 2) this means passing on a simulation object which, iteration by iteration, will grow in size. Using save_cycles in the simulation options only leads to more errors which I will not address here. all in all, I can

Ideally, I would like to be able to activate output_variable in the IDAKLU solver, while retaining the capability of calling and using:

• starting_solution=last_state (in .solve step)
• simulation.solution.cycles[-1]['XXX'].entries (after .solve )
• save_cycles = XXX (in simulation step)

Am I getting this wrong? I.e., recongizning that it might not work just as if I swapped CASADI, are there other workarounds possibile? Let me know :) Thanks for the help

Alessio Lodge, TNO - Imperial

Steps to Reproduce

Code is above. Use any experiment and run with this:

model = pybamm.lithium_ion.DFN()
solver_CasadiSolverSafe = pybamm.CasadiSolver(dt_max = 100,mode = "safe", rtol=1e-7 ,atol=1e-7)
solver_IDAKLU1 = pybamm.IDAKLUSolver(rtol=0.00001, atol=1e-07, root_method='casadi', root_tol=1e-07, extrap_tol=None)
solver_IDAKLU2 = pybamm.IDAKLUSolver(rtol=0.00001, atol=1e-07, root_method='casadi', root_tol=1e-07, extrap_tol=None, output_variables=['Time [s]', 'Current [A]', "Voltage [V]" ], options=None)
solver_IDAKLU3 = pybamm.IDAKLUSolver(rtol=0.00001, atol=1e-07, root_method='casadi', root_tol=1e-07, extrap_tol=None, output_variables=['Time [s]', 'Current [A]', "Voltage [V]", 'Throughput capacity [A.h]', 'Discharge capacity [A.h]'], options= {"dt_max": 1000.0})

FOR CASADI OR IDAKLU WITHOUT OUTPUT VARS

simulation_IDA3_checkup.solution.cycles[-1]['Throughput capacity [A.h]'].entries

FOR IDAKLU WITH OUTPUT VARS 1)

capacity_full_checkup = simulation_IDA3_checkup.solution.get_data_dict()['Throughput capacity [A.h]']
cycles = simulation_IDA3_checkup.solution.get_data_dict()['Cycle']
capacity_discharge_section= [d for d, c in zip(capacity_full_checkup, cycles) if c == 2]   

#I find this work around inconvenient and would prefer that what works with CASADI and IDAKLU (without output_Variables selected) works also with #output_variables selected

FOR IDAKLU WITH OUTPUT VARS 2)

starting_solution = simulation_IDA3_ageing_calendar.solution #

Relevant log output

for 1) (when I make the error happen, without work around)

  File "c:\Users\XXX\Desktop\Journal pape\Code\Python_ageing_code_working_in_IDAKLU.py", line 81, in get_EOL_cap
    print(simulation_IDA3_checkup.solution.cycles[-1]['Throughput capacity [A.h]'].entries - min(simulation_IDA3_checkup.solution.cycles[-1]['Throughput capacity [A.h]'].entries)) #DOESNT WORK WITH IDAKLU, WHY???                                                                                                                                                                             WITH IDAKLU, WHY???
          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\lodgeaa\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\solvers\solution.py", line 710, in __getitem__
    self.update(key)
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\solvers\solution.py", line 590, in update
    self._update_variable(variable)
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\solvers\solution.py", line 608, in _update_variable
    raise KeyError(
KeyError: "Cannot process variable 'Throughput capacity [A.h]' as it was not part of the solve. Please re-run the solve with `output_variables` set to include this variable."


for 2) (when I make the error happen, without work around)

Traceback (most recent call last):
  File "c:\Users\XXX\Desktop\Journal pape\Code\python_ageing_code_working_in_IDAKLU_with_output_variables", line 173, in <module>
    simulation_IDA3_checkup.solve(starting_solution=last_state)
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\simulation.py", line 864, in solve
    self._solution = self._solution + cycle_solution
                     ~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\solvers\solution.py", line 994, in __add__
    v: self._variables[v]._update(other._variables[v], new_sol)
       ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\solvers\processed_variable_computed.py", line 475, in _update
    new_var = self.__class__(bv, bvc, bvd, new_sol)
              ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\solvers\processed_variable_computed.py", line 88, in __init__
    self.initialise_0D()
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\pybamm\solvers\processed_variable_computed.py", line 202, in initialise_0D
    self._xr_data_array = xr.DataArray(entries, coords=[("t", self.t_pts)])
                          ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\xarray\core\dataarray.py", line 456, in __init__
    coords, dims = _infer_coords_and_dims(data.shape, coords, dims)
                   ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\xarray\core\dataarray.py", line 195, in _infer_coords_and_dims
    _check_coords_dims(shape, new_coords, dims_tuple)
  File "C:\Users\XXX\AppData\Local\Programs\Python\Python312\Lib\site-packages\xarray\core\dataarray.py", line 129, in _check_coords_dims
    raise ValueError(
ValueError: conflicting sizes for dimension 't': length 143 on the data but length 144 on coordinate 't'

[pipliggins]

Hi @aaltno, thanks for raising these issues! I'm going to try and summarise some of your points, let me know if I've misunderstood anything.

1. There's a bug which means you can't access variables from individual cycles when using the output_variables option with the IDAKLU solver. MWE:

import pybamm

experiment = pybamm.Experiment(
    [
        "Discharge at C/10 for 10 hours or until 3.3 V",
        "Rest for 1 hour",
    ]
)

model = pybamm.lithium_ion.DFN()
solver_IDAKLU3 = pybamm.IDAKLUSolver(
    output_variables=[
        "Time [s]",
        "Current [A]",
        "Voltage [V]",
        "Throughput capacity [A.h]",
        "Discharge capacity [A.h]",
    ],
)

sim = pybamm.Simulation(model, solver=solver_IDAKLU3, experiment=experiment)
sol = sim.solve()
var = sol.cycles[-1]["Throughput capacity [A.h]"]

raises a key error:

KeyError: "Cannot process variable 'Throughput capacity [A.h]' as it was not part of the solve. Please re-run the solve with `output_variables` set to include this variable."

2. When combining the output_variables option with starting_solution in solving an experiment simulation there's a bug that doesn't allow you to use just the last state as the starting solution. MWE:

model = pybamm.lithium_ion.DFN()

experiment = pybamm.Experiment(
    [
        "Discharge at C/2 for 10 minutes",
        "Rest for 10 minutes",
    ]
)

solver_IDAKLU3 = pybamm.IDAKLUSolver(
    output_variables=[
        "Time [s]",
        "Current [A]",
        "Voltage [V]",
        "Throughput capacity [A.h]",
        "Discharge capacity [A.h]",
    ],
)

sim = pybamm.Simulation(model, experiment=experiment, solver=solver_IDAKLU3)
solution = sim.solve()

sim2 = pybamm.Simulation(model, experiment=experiment, solver=solver_IDAKLU3)

# This works
# new_solution = sim2.solve(starting_solution=solution)

# This doesn't
new_solution = sim2.solve(starting_solution=solution.last_state)

Relevant error:

  File "/Users/pipliggins/Documents/repos/pybamm-local/venv/lib/python3.11/site-packages/xarray/core/dataarray.py", line 140, in _check_coords_dims
    raise ValueError(
ValueError: conflicting sizes for dimension 't': length 90 on the data but length 91 on coordinate 't'

3. You're finding errors with the save_cycles option not specified here.

[pipliggins]

(1) Is fairly easily fixed, I'm going to link a draft PR where that bug should be solved if you'd like to try it out and check it does what you're expecting.

I'll take a look at (2); @aaltno could you provide a minimum working example and traceback for the errors you allude to for (3)? You can use the code in the PR first to see if some of your issues are fixed first.