Problem when differentiating FunctionParameter with respect to time

[mlopber]

Hello everyone,

I am building a PyBaMM model and I need to differentiate a pybamm.FunctionParameter whith respect to time. I have the following code:

import pybamm
from pybamm import FunctionParameter
from utils_functions import nu_n_exp

model = pybamm.lithium_ion.DFN()
param = pybamm.ParameterValues('Chen2020')

time = model.variables["Time [s]"]
sto_n = model.variables["Negative electrode stoichiometry"]


def nu_n(sto) -> FunctionParameter:
    return FunctionParameter("Negative electrode expansion function [m]",
                             {"Negative electrode stoichiometry": sto})


param.update({"Negative electrode expansion function [m]": nu_n_exp},
             check_already_exists=False)

dnu_n_dt_fun = nu_n(sto=sto_n).diff(variable=time)
# dnu_n_dt_fun = nu_n(sto=sto_n).diff(variable=sto_n)

model.variables.update({
    "Negative electrode expansion rate [m/s]": dnu_n_dt_fun,
})


# experiment
experiment = pybamm.Experiment([("Discharge at 1C for 1 hour")])

# simulaton instance
sim = pybamm.Simulation(model, experiment=experiment, parameter_values=param,
                        solver=pybamm.IDAKLUSolver())

# Solve the model
sim.solve(initial_soc=1)

When I try to run a customized DFN model with this pybamm.FunctionParameter derivative I get the following error:

""" Default behaviour for differentiation, overriden by Binary and Unary Operators """ raise NotImplementedError

Nevertheless, if I differentiate the nu_n variable with respect to "sto" (stoichiometry) I have no problem at all. Do you guys know a way to circumvent this issue?

Kind regards

[martinjrobins]

Looks like pybamm can't take the derivative of one of the nodes in the Negative electrode stoichiometry expression:

x-average(0.75 * r-average(Negative particle concentration [mol.m-3]))

The x-average is the problem here, diff doesn't seem to be implemented for any of the averages. The problem doesn't occur if you take the derivative wrt sto_n because pybamm can work that out without taking the derivative of the average (since d sto /d sto = 1)

You could circumvent this if you manually work out the derivative of nu_n and give this to a pybamm.FunctionParameter. Mind you, if the result contains the time derivative of one of the state variables (e.g. concentration) then you might not be able to solve the model at all, since we assume no time derivatives of the state variables in the RHS of the equations. If this is the case then you will need to re-structure your model so you are back to the form dy / dt = f(y, t, p)