Skip to content

Commit

Permalink
Merge branch 'develop' into issue-4183-remove-autograd
Browse files Browse the repository at this point in the history
  • Loading branch information
brosaplanella authored Jun 21, 2024
2 parents 84ee101 + 8ba4791 commit d38117b
Show file tree
Hide file tree
Showing 22 changed files with 402 additions and 59 deletions.
9 changes: 9 additions & 0 deletions .all-contributorsrc
Original file line number Diff line number Diff line change
Expand Up @@ -922,6 +922,15 @@
"code",
"infra"
]
},
{
"login": "smitasahu2",
"name": "Smita Sahu",
"avatar_url": "https://avatars.githubusercontent.com/u/57876346?v=4",
"profile": "https://github.com/smitasahu2",
"contributions": [
"code"
]
}
],
"contributorsPerLine": 7,
Expand Down
6 changes: 4 additions & 2 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@

- Added new parameters `"f{pref]Initial inner SEI on cracks thickness [m]"` and `"f{pref]Initial outer SEI on cracks thickness [m]"`, instead of hardcoding these to `L_inner_0 / 10000` and `L_outer_0 / 10000`. ([#4168](https://github.com/pybamm-team/PyBaMM/pull/4168))
- Added `pybamm.DataLoader` class to fetch data files from [pybamm-data](https://github.com/pybamm-team/pybamm-data/releases/tag/v1.0.0) and store it under local cache. ([#4098](https://github.com/pybamm-team/PyBaMM/pull/4098))
- Transport efficiency submodel has new options from the literature relating to different tortuosity factor models and also a new option called "tortuosity factor" for specifying the value or function directly as parameters ([#3437](https://github.com/pybamm-team/PyBaMM/pull/3437))
- Added `time` as an option for `Experiment.termination`. Now allows solving up to a user-specified time while also allowing different cycles and steps in an experiment to be handled normally. ([#4073](https://github.com/pybamm-team/PyBaMM/pull/4073))
- Added `plot_thermal_components` to plot the contributions to the total heat generation in a battery ([#4021](https://github.com/pybamm-team/PyBaMM/pull/4021))
- Added functions for normal probability density function (`pybamm.normal_pdf`) and cumulative distribution function (`pybamm.normal_cdf`) ([#3999](https://github.com/pybamm-team/PyBaMM/pull/3999))
- "Basic" models are now compatible with experiments ([#3995](https://github.com/pybamm-team/PyBaMM/pull/3995))
Expand All @@ -20,10 +20,12 @@
- Add support for BPX version 0.4.0 which allows for blended electrodes and user-defined parameters in BPX([#3414](https://github.com/pybamm-team/PyBaMM/pull/3414))
- Added `by_submodel` feature in `print_parameter_info` method to allow users to print parameters and types of submodels in a tabular and readable format ([#3628](https://github.com/pybamm-team/PyBaMM/pull/3628))
- Added `WyciskOpenCircuitPotential` for differential capacity hysteresis state open-circuit potential submodel ([#3593](https://github.com/pybamm-team/PyBaMM/pull/3593))
- Added `time` as an option for `Experiment.termination`. Now allows solving up to a user-specified time while also allowing different cycles and steps in an experiment to be handled normally. ([#4073](https://github.com/pybamm-team/PyBaMM/pull/4073))
- Transport efficiency submodel has new options from the literature relating to different tortuosity factor models and also a new option called "tortuosity factor" for specifying the value or function directly as parameters ([#3437](https://github.com/pybamm-team/PyBaMM/pull/3437))
- Heat of mixing source term can now be included into thermal models ([#2837](https://github.com/pybamm-team/PyBaMM/pull/2837))

## Bug Fixes

- Fixed bug where passing deprecated `electrode diffusivity` parameter resulted in a breaking change and/or the corresponding diffusivity parameter not updating. Improved the deprecated translation around BPX. ([#4176](https://github.com/pybamm-team/PyBaMM/pull/4176))
- Fixed a bug where a factor of electrode surface area to volume ratio is missing in the rhs of the LeadingOrderDifferential conductivity model ([#4139](https://github.com/pybamm-team/PyBaMM/pull/4139))
- Fixes the breaking changes caused by [#3624](https://github.com/pybamm-team/PyBaMM/pull/3624), specifically enables the deprecated parameter `electrode diffusivity` to be used by `ParameterValues.update({name:value})` and `Solver.solve(inputs={name:value})`. Fixes parameter translation from old name to new name, with corrected tests. ([#4072](https://github.com/pybamm-team/PyBaMM/pull/4072)
- Set the `remove_independent_variables_from_rhs` to `False` by default, and moved the option from `Discretisation.process_model` to `Discretisation.__init__`. This fixes a bug related to the discharge capacity, but may make the simulation slower in some cases. To set the option to `True`, use `Simulation(..., discretisation_kwargs={"remove_independent_variables_from_rhs": True})`. ([#4020](https://github.com/pybamm-team/PyBaMM/pull/4020))
Expand Down
2 changes: 1 addition & 1 deletion README.md
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@
[![code style](https://img.shields.io/endpoint?url=https://raw.githubusercontent.com/charliermarsh/ruff/main/assets/badge/v2.json)](https://github.com/astral-sh/ruff)

<!-- ALL-CONTRIBUTORS-BADGE:START - Do not remove or modify this section -->
[![All Contributors](https://img.shields.io/badge/all_contributors-86-orange.svg)](#-contributors)
[![All Contributors](https://img.shields.io/badge/all_contributors-87-orange.svg)](#-contributors)
<!-- ALL-CONTRIBUTORS-BADGE:END -->

</div>
Expand Down
1 change: 1 addition & 0 deletions all_contributors.md
Original file line number Diff line number Diff line change
Expand Up @@ -116,6 +116,7 @@ Thanks goes to these wonderful people ([emoji key](https://allcontributors.org/d
<tr>
<td align="center" valign="top" width="14.28%"><a href="https://github.com/ikorotkin"><img src="https://avatars.githubusercontent.com/u/29599800?v=4?s=100" width="100px;" alt="Ivan Korotkin"/><br /><sub><b>Ivan Korotkin</b></sub></a><br /><a href="https://github.com/pybamm-team/PyBaMM/commits?author=ikorotkin" title="Code">💻</a></td>
<td align="center" valign="top" width="14.28%"><a href="https://github.com/santacodes"><img src="https://avatars.githubusercontent.com/u/52504160?v=4?s=100" width="100px;" alt="Santhosh"/><br /><sub><b>Santhosh</b></sub></a><br /><a href="https://github.com/pybamm-team/PyBaMM/commits?author=santacodes" title="Code">💻</a> <a href="#infra-santacodes" title="Infrastructure (Hosting, Build-Tools, etc)">🚇</a></td>
<td align="center" valign="top" width="14.28%"><a href="https://github.com/smitasahu2"><img src="https://avatars.githubusercontent.com/u/57876346?v=4?s=100" width="100px;" alt="Smita Sahu"/><br /><sub><b>Smita Sahu</b></sub></a><br /><a href="https://github.com/pybamm-team/PyBaMM/commits?author=smitasahu2" title="Code">💻</a></td>
</tr>
</tbody>
</table>
Expand Down
112 changes: 112 additions & 0 deletions examples/scripts/compare_lithium_ion_heat_of_mixing.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,112 @@
#
# Compare SPMe model with and without heat of mixing
#
import pybamm
import numpy as np
import matplotlib.pyplot as plt

pybamm.set_logging_level("INFO")

# load models
models = [
pybamm.lithium_ion.SPMe(
{"heat of mixing": "true", "thermal": "lumped"}, name="SPMe with heat of mixing"
),
pybamm.lithium_ion.SPMe({"thermal": "lumped"}, name="SPMe without heat of mixing"),
]

# set parametrisation (Ecker et al., 2015)
parameter_values = pybamm.ParameterValues("Ecker2015")

# set mesh
# (increased number of points in particles to avoid oscillations for high C-rates)
var_pts = {"x_n": 16, "x_s": 8, "x_p": 16, "r_n": 128, "r_p": 128}

# set the constant current discharge C-rate
C_rate = 5

# set the simulation time and increase the number of points
# for better integration in time
t_eval = np.linspace(0, 720, 360)

# set up an extra plot with the heat of mixing vs time in each electrode and
# the integrated heat of mixing vs time in each electrode to compare with
# Figure 6(a) from Richardson et al. (2021)
fig, axs = plt.subplots(2, len(models), figsize=(12, 7))

# extract some of the parameters
L_n = parameter_values["Negative electrode thickness [m]"]
L_p = parameter_values["Positive electrode thickness [m]"]
A = parameter_values["Electrode width [m]"] * parameter_values["Electrode height [m]"]

# create and run simulations
sims = []
for m, model in enumerate(models):
sim = pybamm.Simulation(
model, parameter_values=parameter_values, var_pts=var_pts, C_rate=C_rate
)
sim.solve(t_eval)
sims.append(sim)

# extract solution for plotting
sol = sim.solution

# extract variables from the solution
time = sol["Time [h]"].entries
Q_mix = sol["Heat of mixing [W.m-3]"].entries

# heat of mixing in negative and positive electrodes multiplied by the electrode
# width, represents the integral of heat of mixing term across each of the
# electrodes (W.m-2)
Q_mix_n = Q_mix[0, :] * L_n
Q_mix_p = Q_mix[-1, :] * L_p

# heat of mixing integrals (J.m-2)
Q_mix_n_int = 0.0
Q_mix_p_int = 0.0

# data for plotting
Q_mix_n_plt = []
Q_mix_p_plt = []

# performs integration in time
for i, t in enumerate(time[1:]):
dt = (t - time[i]) * 3600 # seconds
Q_mix_n_avg = (Q_mix_n[i] + Q_mix_n[i + 1]) * 0.5
Q_mix_p_avg = (Q_mix_p[i] + Q_mix_p[i + 1]) * 0.5
# convert J to kJ and divide the integral by the electrode area A to compare
# with Figure 6(a) from Richardson et al. (2021)
Q_mix_n_int += Q_mix_n_avg * dt / 1000 / A
Q_mix_p_int += Q_mix_p_avg * dt / 1000 / A
Q_mix_n_plt.append(Q_mix_n_int)
Q_mix_p_plt.append(Q_mix_p_int)

# plots heat of mixing in each electrode vs time in minutes
axs[0, m].plot(time * 60, Q_mix_n, ls="-", label="Negative electrode")
axs[0, m].plot(time * 60, Q_mix_p, ls="--", label="Positive electrode")
axs[0, m].set_title(f"{model.name}")
axs[0, m].set_xlabel("Time [min]")
axs[0, m].set_ylabel("Heat of mixing [W.m-2]")
axs[0, m].grid(True)
axs[0, m].legend()

# plots integrated heat of mixing in each electrode vs time in minutes
axs[1, m].plot(time[1:] * 60, Q_mix_n_plt, ls="-", label="Negative electrode")
axs[1, m].plot(time[1:] * 60, Q_mix_p_plt, ls="--", label="Positive electrode")
axs[1, m].set_xlabel("Time [min]")
axs[1, m].set_ylabel("Integrated heat of mixing [kJ.m-2]")
axs[1, m].grid(True)
axs[1, m].legend()

# plot
pybamm.dynamic_plot(
sims,
output_variables=[
"X-averaged cell temperature [K]",
"X-averaged heat of mixing [W.m-3]",
"X-averaged total heating [W.m-3]",
"Heat of mixing [W.m-3]",
"Voltage [V]",
"Current [A]",
],
)
11 changes: 11 additions & 0 deletions pybamm/CITATIONS.bib
Original file line number Diff line number Diff line change
Expand Up @@ -418,6 +418,17 @@ @article{Richardson2020
doi = {10.1016/j.electacta.2020.135862},
}

@article{Richardson2021,
title = {Heat Generation and a Conservation Law for Chemical Energy in Li-ion Batteries},
author = {Richardson, Giles and Korotkin, Ivan},
journal = {Electrochimica Acta},
volume = {392},
pages = {138909},
year = {2021},
publisher = {Elsevier},
doi = {10.1016/j.electacta.2021.138909},
}

@article{Sripad2020,
title={Kinetics of lithium electrodeposition and stripping},
author={Sripad, Shashank and Korff, Daniel and DeCaluwe, Steven C and Viswanathan, Venkatasubramanian},
Expand Down
11 changes: 10 additions & 1 deletion pybamm/models/full_battery_models/base_battery_model.py
Original file line number Diff line number Diff line change
Expand Up @@ -238,6 +238,7 @@ def __init__(self, extra_options):
"integrated",
],
"exchange-current density": ["single", "current sigmoid"],
"heat of mixing": ["false", "true"],
"hydrolysis": ["false", "true"],
"intercalation kinetics": [
"symmetric Butler-Volmer",
Expand Down Expand Up @@ -512,6 +513,11 @@ def __init__(self, extra_options):

# Options not yet compatible with particle-size distributions
if options["particle size"] == "distribution":
if options["heat of mixing"] != "false":
raise NotImplementedError(
"Heat of mixing submodels do not yet support particle-size "
"distributions."
)
if options["lithium plating"] != "none":
raise NotImplementedError(
"Lithium plating submodels do not yet support particle-size "
Expand Down Expand Up @@ -1240,7 +1246,10 @@ def set_thermal_submodel(self):
if self.options["dimensionality"] == 0:
thermal_submodel = pybamm.thermal.pouch_cell.OneDimensionalX

self.submodels["thermal"] = thermal_submodel(self.param, self.options)
x_average = getattr(self, "x_average", False)
self.submodels["thermal"] = thermal_submodel(
self.param, self.options, x_average
)

def set_current_collector_submodel(self):
if self.options["current collector"] in ["uniform"]:
Expand Down
94 changes: 92 additions & 2 deletions pybamm/models/submodels/thermal/base_thermal.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
# Base class for thermal effects
#
import pybamm
import numpy as np


class BaseThermal(pybamm.BaseSubModel):
Expand All @@ -16,8 +17,12 @@ class BaseThermal(pybamm.BaseSubModel):
A dictionary of options to be passed to the model.
"""

def __init__(self, param, options=None):
def __init__(self, param, options=None, x_average=False):
super().__init__(param, options=options)
self.x_average = x_average

if self.options["heat of mixing"] == "true":
pybamm.citations.register("Richardson2021")

def _get_standard_fundamental_variables(self, T_dict):
"""
Expand Down Expand Up @@ -175,15 +180,20 @@ def _get_standard_coupled_variables(self, variables):
Q_rev_n, pybamm.FullBroadcast(0, "separator", "current collector"), Q_rev_p
)

# Heat of mixing
Q_mix_s_n, Q_mix_s_s, Q_mix_s_p = self._heat_of_mixing(variables)
Q_mix = pybamm.concatenation(Q_mix_s_n, Q_mix_s_s, Q_mix_s_p)

# Total heating
Q = Q_ohm + Q_rxn + Q_rev
Q = Q_ohm + Q_rxn + Q_rev + Q_mix

# Compute the X-average over the entire cell, including current collectors
# Note: this can still be a function of y and z for higher-dimensional pouch
# cell models
Q_ohm_av = self._x_average(Q_ohm, Q_ohm_s_cn, Q_ohm_s_cp)
Q_rxn_av = self._x_average(Q_rxn, 0, 0)
Q_rev_av = self._x_average(Q_rev, 0, 0)
Q_mix_av = self._x_average(Q_mix, 0, 0)
Q_av = self._x_average(Q, Q_ohm_s_cn, Q_ohm_s_cp)

# Compute the integrated heat source per unit simulated electrode-pair area
Expand All @@ -192,11 +202,14 @@ def _get_standard_coupled_variables(self, variables):
Q_ohm_Wm2 = Q_ohm_av * param.L
Q_rxn_Wm2 = Q_rxn_av * param.L
Q_rev_Wm2 = Q_rev_av * param.L
Q_mix_Wm2 = Q_mix_av * param.L
Q_Wm2 = Q_av * param.L

# Now average over the electrode height and width
Q_ohm_Wm2_av = self._yz_average(Q_ohm_Wm2)
Q_rxn_Wm2_av = self._yz_average(Q_rxn_Wm2)
Q_rev_Wm2_av = self._yz_average(Q_rev_Wm2)
Q_mix_Wm2_av = self._yz_average(Q_mix_Wm2)
Q_Wm2_av = self._yz_average(Q_Wm2)

# Compute total heat source terms (in W) over the *entire cell volume*, not
Expand All @@ -208,6 +221,7 @@ def _get_standard_coupled_variables(self, variables):
Q_ohm_W = Q_ohm_Wm2_av * n_elec * A
Q_rxn_W = Q_rxn_Wm2_av * n_elec * A
Q_rev_W = Q_rev_Wm2_av * n_elec * A
Q_mix_W = Q_mix_Wm2_av * n_elec * A
Q_W = Q_Wm2_av * n_elec * A

# Compute volume-averaged heat source terms over the *entire cell volume*, not
Expand All @@ -216,6 +230,7 @@ def _get_standard_coupled_variables(self, variables):
Q_ohm_vol_av = Q_ohm_W / V
Q_rxn_vol_av = Q_rxn_W / V
Q_rev_vol_av = Q_rev_W / V
Q_mix_vol_av = Q_mix_W / V
Q_vol_av = Q_W / V

# Effective heat capacity
Expand All @@ -228,6 +243,7 @@ def _get_standard_coupled_variables(self, variables):
"Ohmic heating [W.m-3]": Q_ohm,
"X-averaged Ohmic heating [W.m-3]": Q_ohm_av,
"Volume-averaged Ohmic heating [W.m-3]": Q_ohm_vol_av,
"Volume-averaged heat of mixing [W.m-3]": Q_mix_vol_av,
"Ohmic heating per unit electrode-pair area [W.m-2]": Q_ohm_Wm2,
"Ohmic heating [W]": Q_ohm_W,
# Irreversible
Expand All @@ -244,6 +260,12 @@ def _get_standard_coupled_variables(self, variables):
"Volume-averaged reversible heating [W.m-3]": Q_rev_vol_av,
"Reversible heating per unit electrode-pair area " "[W.m-2]": Q_rev_Wm2,
"Reversible heating [W]": Q_rev_W,
# Mixing
"Heat of mixing [W.m-3]": Q_mix,
"X-averaged heat of mixing [W.m-3]": Q_mix_av,
"Volume-averaged heating of mixing [W.m-3]": Q_mix_vol_av,
"Heat of mixing per unit electrode-pair area " "[W.m-2]": Q_mix_Wm2,
"Heat of mixing [W]": Q_mix_W,
# Total
"Total heating [W.m-3]": Q,
"X-averaged total heating [W.m-3]": Q_av,
Expand Down Expand Up @@ -290,6 +312,74 @@ def _current_collector_heating(self, variables):
Q_s_cp = self.param.p.sigma_cc * pybamm.grad_squared(phi_s_cp)
return Q_s_cn, Q_s_cp

def _heat_of_mixing(self, variables):
"""Compute heat of mixing source terms."""
param = self.param

if self.options["heat of mixing"] == "true":
F = pybamm.constants.F.value
pi = np.pi

# Compute heat of mixing in negative electrode
if self.options.electrode_types["negative"] == "planar":
Q_mix_s_n = pybamm.FullBroadcast(
0, ["negative electrode"], "current collector"
)
else:
a_n = variables["Negative electrode surface area to volume ratio [m-1]"]
R_n = variables["Negative particle radius [m]"]
N_n = a_n / (4 * pi * R_n**2)
if self.x_average:
c_n = variables[
"X-averaged negative particle concentration [mol.m-3]"
]
T_n = variables["X-averaged negative electrode temperature [K]"]
else:
c_n = variables["Negative particle concentration [mol.m-3]"]
T_n = variables["Negative electrode temperature [K]"]
T_n_part = pybamm.PrimaryBroadcast(T_n, ["negative particle"])
dc_n_dr2 = pybamm.inner(pybamm.grad(c_n), pybamm.grad(c_n))
D_n = param.n.prim.D(c_n, T_n_part)
dUeq_n = param.n.prim.U(c_n / param.n.prim.c_max, T_n_part).diff(c_n)
integrand_r_n = D_n * dc_n_dr2 * dUeq_n
integration_variable_r_n = [
pybamm.SpatialVariable("r", domain=integrand_r_n.domain)
]
integral_r_n = pybamm.Integral(integrand_r_n, integration_variable_r_n)
Q_mix_s_n = -F * N_n * integral_r_n

# Compute heat of mixing in positive electrode
a_p = variables["Positive electrode surface area to volume ratio [m-1]"]
R_p = variables["Positive particle radius [m]"]
N_p = a_p / (4 * pi * R_p**2)
if self.x_average:
c_p = variables["X-averaged positive particle concentration [mol.m-3]"]
T_p = variables["X-averaged positive electrode temperature [K]"]
else:
c_p = variables["Positive particle concentration [mol.m-3]"]
T_p = variables["Positive electrode temperature [K]"]
T_p_part = pybamm.PrimaryBroadcast(T_p, ["positive particle"])
dc_p_dr2 = pybamm.inner(pybamm.grad(c_p), pybamm.grad(c_p))
D_p = param.p.prim.D(c_p, T_p_part)
dUeq_p = param.p.prim.U(c_p / param.p.prim.c_max, T_p_part).diff(c_p)
integrand_r_p = D_p * dc_p_dr2 * dUeq_p
integration_variable_r_p = [
pybamm.SpatialVariable("r", domain=integrand_r_p.domain)
]
integral_r_p = pybamm.Integral(integrand_r_p, integration_variable_r_p)
Q_mix_s_p = -F * N_p * integral_r_p
Q_mix_s_s = pybamm.FullBroadcast(0, ["separator"], "current collector")
else:
Q_mix_s_n = pybamm.FullBroadcast(
0, ["negative electrode"], "current collector"
)
Q_mix_s_p = pybamm.FullBroadcast(
0, ["positive electrode"], "current collector"
)
Q_mix_s_s = pybamm.FullBroadcast(0, ["separator"], "current collector")

return Q_mix_s_n, Q_mix_s_s, Q_mix_s_p

def _x_average(self, var, var_cn, var_cp):
"""
Computes the X-average over the whole cell (including current collectors)
Expand Down
Loading

0 comments on commit d38117b

Please sign in to comment.