diff --git a/CHANGELOG.md b/CHANGELOG.md index b2f8637442..57e0aad007 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,14 +1,15 @@ # [Unreleased](https://github.com/pybamm-team/PyBaMM/) -- Replaced rounded faraday constant with its exact value in `bpx.py` for better comparison between different tools - ## Features - Added additional user-configurable options to the (`IDAKLUSolver`) and adjusted the default values to improve performance. ([#4282](https://github.com/pybamm-team/PyBaMM/pull/4282)) +- Added the diffusion element to be used in the Thevenin model. ([#4254](https://github.com/pybamm-team/PyBaMM/pull/4254)) ## Optimizations - Improved performance and reliability of DAE consistent initialization. ([#4301](https://github.com/pybamm-team/PyBaMM/pull/4301)) +- Replaced rounded Faraday constant with its exact value in `bpx.py` for better comparison between different tools. ([#4290](https://github.com/pybamm-team/PyBaMM/pull/4290)) + ## Bug Fixes - Fixed bug where IDAKLU solver failed when `output variables` were specified and an event triggered. ([#4300](https://github.com/pybamm-team/PyBaMM/pull/4300)) diff --git a/examples/scripts/run_ecmd.py b/examples/scripts/run_ecmd.py new file mode 100644 index 0000000000..c12fef9cfb --- /dev/null +++ b/examples/scripts/run_ecmd.py @@ -0,0 +1,28 @@ +import pybamm + +pybamm.set_logging_level("INFO") + +model = pybamm.equivalent_circuit.Thevenin(options={"diffusion element": "true"}) +parameter_values = model.default_parameter_values + +parameter_values.update( + {"Diffusion time constant [s]": 580}, check_already_exists=False +) + +experiment = pybamm.Experiment( + [ + ( + "Discharge at C/10 for 10 hours or until 3.3 V", + "Rest for 30 minutes", + "Rest for 2 hours", + "Charge at 100 A until 4.1 V", + "Hold at 4.1 V until 5 A", + "Rest for 30 minutes", + "Rest for 1 hour", + ), + ] +) + +sim = pybamm.Simulation(model, experiment=experiment, parameter_values=parameter_values) +sim.solve() +sim.plot() diff --git a/pybamm/CITATIONS.bib b/pybamm/CITATIONS.bib index 33ac6055d3..8fb9c6dc98 100644 --- a/pybamm/CITATIONS.bib +++ b/pybamm/CITATIONS.bib @@ -794,3 +794,13 @@ @article{Wycisk2022 author = {Dominik Wycisk and Marc Oldenburger and Marc Gerry Stoye and Toni Mrkonjic and Arnulf Latz}, keywords = {Lithium-ion battery, Voltage hysteresis, Plett-model, Silicon–graphite anode}, } + +@article{Fan2022, + title={Data-driven identification of lithium-ion batteries: A nonlinear equivalent circuit model with diffusion dynamics}, + author={Fan, Chuanxin and O’Regan, Kieran and Li, Liuying and Higgins, Matthew D and Kendrick, Emma and Widanage, Widanalage D}, + journal={Applied Energy}, + volume={321}, + pages={119336}, + year={2022}, + publisher={Elsevier} +} diff --git a/pybamm/models/full_battery_models/equivalent_circuit/thevenin.py b/pybamm/models/full_battery_models/equivalent_circuit/thevenin.py index 9aaf747a4c..48c0c1d1ac 100644 --- a/pybamm/models/full_battery_models/equivalent_circuit/thevenin.py +++ b/pybamm/models/full_battery_models/equivalent_circuit/thevenin.py @@ -30,6 +30,9 @@ class Thevenin(pybamm.BaseModel): throughput capacity in addition to discharge capacity. Must be one of "true" or "false". "false" is the default, since calculating discharge energy can be computationally expensive for simple models like SPM. + * "diffusion element" : str + Whether to include the diffusion element to the model. Must be one of + "true" or "false". "false" is the default. * "operating mode" : str Sets the operating mode for the model. This determines how the current is set. Can be: @@ -73,6 +76,7 @@ def __init__( def set_options(self, extra_options=None): possible_options = { "calculate discharge energy": ["false", "true"], + "diffusion element": ["false", "true"], "operating mode": OperatingModes("current"), "number of rc elements": NaturalNumberOption(1), } @@ -165,6 +169,18 @@ def set_rc_submodels(self): ) self.element_counter += 1 + def set_diffusion_submodel(self): + if self.options["diffusion element"] == "false": + self.submodels["Diffusion"] = ( + pybamm.equivalent_circuit_elements.NoDiffusion(self.param, self.options) + ) + elif self.options["diffusion element"] == "true": + self.submodels["Diffusion"] = ( + pybamm.equivalent_circuit_elements.DiffusionElement( + self.param, self.options + ) + ) + def set_thermal_submodel(self): self.submodels["Thermal"] = pybamm.equivalent_circuit_elements.ThermalSubModel( self.param, self.options @@ -180,6 +196,7 @@ def set_submodels(self, build): self.set_ocv_submodel() self.set_resistor_submodel() self.set_rc_submodels() + self.set_diffusion_submodel() self.set_thermal_submodel() self.set_voltage_submodel() @@ -227,3 +244,25 @@ def default_quick_plot_variables(self): "Irreversible heat generation [W]", ], ] + + @property + def default_var_pts(self): + x = pybamm.SpatialVariable( + "x ECMD", domain=["ECMD particle"], coord_sys="Cartesian" + ) + return {x: 20} + + @property + def default_geometry(self): + x = pybamm.SpatialVariable( + "x ECMD", domain=["ECMD particle"], coord_sys="Cartesian" + ) + return {"ECMD particle": {x: {"min": 0, "max": 1}}} + + @property + def default_submesh_types(self): + return {"ECMD particle": pybamm.Uniform1DSubMesh} + + @property + def default_spatial_methods(self): + return {"ECMD particle": pybamm.FiniteVolume()} diff --git a/pybamm/models/submodels/equivalent_circuit_elements/__init__.py b/pybamm/models/submodels/equivalent_circuit_elements/__init__.py index 4ff6ee7c62..16d57c691c 100644 --- a/pybamm/models/submodels/equivalent_circuit_elements/__init__.py +++ b/pybamm/models/submodels/equivalent_circuit_elements/__init__.py @@ -3,6 +3,7 @@ from .rc_element import RCElement from .thermal import ThermalSubModel from .voltage_model import VoltageModel +from .diffusion_element import NoDiffusion, DiffusionElement __all__ = ['ocv_element', 'rc_element', 'resistor_element', 'thermal', 'voltage_model'] diff --git a/pybamm/models/submodels/equivalent_circuit_elements/diffusion_element.py b/pybamm/models/submodels/equivalent_circuit_elements/diffusion_element.py new file mode 100644 index 0000000000..c7f1b4bcd5 --- /dev/null +++ b/pybamm/models/submodels/equivalent_circuit_elements/diffusion_element.py @@ -0,0 +1,110 @@ +import pybamm + + +class NoDiffusion(pybamm.BaseSubModel): + """ + Without Diffusion element for + equivalent circuits. + + Parameters + ---------- + param : parameter class + The parameters to use for this submodel + options : dict, optional + A dictionary of options to be passed to the model. + """ + + def __init__(self, param, options=None): + super().__init__(param) + self.model_options = options + + def get_coupled_variables(self, variables): + z = pybamm.PrimaryBroadcast(variables["SoC"], "ECMD particle") + x = pybamm.SpatialVariable( + "x ECMD", domain=["ECMD particle"], coord_sys="Cartesian" + ) + z_surf = pybamm.surf(z) + eta_diffusion = pybamm.Scalar(0) + + variables.update( + { + "Distributed SoC": z, + "x ECMD": x, + "Diffusion overpotential [V]": eta_diffusion, + "Surface SoC": z_surf, + } + ) + + return variables + + +class DiffusionElement(pybamm.BaseSubModel): + """ + With Diffusion element for + equivalent circuits. + + Parameters + ---------- + param : parameter class + The parameters to use for this submodel + options : dict, optional + A dictionary of options to be passed to the model. + """ + + def __init__(self, param, options=None): + super().__init__(param) + pybamm.citations.register("Fan2022") + self.model_options = options + + def get_fundamental_variables(self): + z = pybamm.Variable("Distributed SoC", domain="ECMD particle") + x = pybamm.SpatialVariable( + "x ECMD", domain=["ECMD particle"], coord_sys="Cartesian" + ) + variables = { + "Distributed SoC": z, + "x ECMD": x, + } + return variables + + def get_coupled_variables(self, variables): + z = variables["Distributed SoC"] + soc = variables["SoC"] + z_surf = pybamm.surf(z) + eta_diffusion = -(self.param.ocv(z_surf) - self.param.ocv(soc)) + + variables.update( + { + "Diffusion overpotential [V]": eta_diffusion, + "Surface SoC": z_surf, + } + ) + + return variables + + def set_rhs(self, variables): + cell_capacity = self.param.cell_capacity + current = variables["Current [A]"] + z = variables["Distributed SoC"] + + # governing equations + dzdt = pybamm.div(pybamm.grad(z)) / self.param.tau_D + self.rhs = {z: dzdt} + + # boundary conditions + lbc = pybamm.Scalar(0) + rbc = -self.param.tau_D * current / (cell_capacity * 3600) + self.boundary_conditions = { + z: {"left": (lbc, "Neumann"), "right": (rbc, "Neumann")} + } + + def set_initial_conditions(self, variables): + z = variables["Distributed SoC"] + self.initial_conditions = {z: self.param.initial_soc} + + def set_events(self, variables): + z_surf = variables["Surface SoC"] + self.events += [ + pybamm.Event("Minimum surface SoC", z_surf), + pybamm.Event("Maximum surface SoC", 1 - z_surf), + ] diff --git a/pybamm/models/submodels/equivalent_circuit_elements/ocv_element.py b/pybamm/models/submodels/equivalent_circuit_elements/ocv_element.py index 4c570b3cf4..9d1adf3d57 100644 --- a/pybamm/models/submodels/equivalent_circuit_elements/ocv_element.py +++ b/pybamm/models/submodels/equivalent_circuit_elements/ocv_element.py @@ -56,7 +56,7 @@ def set_initial_conditions(self, variables): def set_events(self, variables): soc = variables["SoC"] - self.events = [ + self.events += [ pybamm.Event("Minimum SoC", soc), pybamm.Event("Maximum SoC", 1 - soc), ] diff --git a/pybamm/models/submodels/equivalent_circuit_elements/voltage_model.py b/pybamm/models/submodels/equivalent_circuit_elements/voltage_model.py index d976cd0747..380902fca5 100644 --- a/pybamm/models/submodels/equivalent_circuit_elements/voltage_model.py +++ b/pybamm/models/submodels/equivalent_circuit_elements/voltage_model.py @@ -30,7 +30,9 @@ def get_coupled_variables(self, variables): for i in range(number_of_elements): overpotential += variables[f"Element-{i} overpotential [V]"] - voltage = ocv + overpotential + diffusion_overpotential = variables["Diffusion overpotential [V]"] + + voltage = ocv + overpotential + diffusion_overpotential # Power and Resistance current = variables["Current [A]"] diff --git a/pybamm/parameters/ecm_parameters.py b/pybamm/parameters/ecm_parameters.py index 79bb42e318..4d5855c2fd 100644 --- a/pybamm/parameters/ecm_parameters.py +++ b/pybamm/parameters/ecm_parameters.py @@ -4,6 +4,7 @@ class EcmParameters: def __init__(self): self.cell_capacity = pybamm.Parameter("Cell capacity [A.h]") + self.tau_D = pybamm.Parameter("Diffusion time constant [s]") self._set_current_parameters() self._set_voltage_parameters() diff --git a/tests/integration/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py b/tests/integration/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py index e01cf1a7d7..c8c68ea7f0 100644 --- a/tests/integration/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py +++ b/tests/integration/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py @@ -7,3 +7,15 @@ def test_basic_processing(self): model = pybamm.equivalent_circuit.Thevenin() modeltest = tests.StandardModelTest(model) modeltest.test_all() + + def test_diffusion(self): + model = pybamm.equivalent_circuit.Thevenin( + options={"diffusion element": "true"} + ) + parameter_values = model.default_parameter_values + + parameter_values.update( + {"Diffusion time constant [s]": 580}, check_already_exists=False + ) + modeltest = tests.StandardModelTest(model, parameter_values=parameter_values) + modeltest.test_all() diff --git a/tests/unit/test_citations.py b/tests/unit/test_citations.py index ba216e62ff..978569d864 100644 --- a/tests/unit/test_citations.py +++ b/tests/unit/test_citations.py @@ -345,6 +345,18 @@ def test_msmr(self): self.assertIn("Verbrugge2017", citations._papers_to_cite) self.assertIn("Verbrugge2017", citations._citation_tags.keys()) + def test_thevenin(self): + citations = pybamm.citations + + citations._reset() + pybamm.equivalent_circuit.Thevenin() + self.assertNotIn("Fan2022", citations._papers_to_cite) + self.assertNotIn("Fan2022", citations._citation_tags.keys()) + + pybamm.equivalent_circuit.Thevenin(options={"diffusion element": "true"}) + self.assertIn("Fan2022", citations._papers_to_cite) + self.assertIn("Fan2022", citations._citation_tags.keys()) + def test_parameter_citations(self): citations = pybamm.citations diff --git a/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py b/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py index f6abf592f8..39bfaf1145 100644 --- a/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py +++ b/tests/unit/test_models/test_full_battery_models/test_equivalent_circuit/test_thevenin.py @@ -11,6 +11,28 @@ def test_standard_model(self): model = pybamm.equivalent_circuit.Thevenin() model.check_well_posedness() + def test_default_properties(self): + model = pybamm.equivalent_circuit.Thevenin() + x = model.variables["x ECMD"] + + # test var_pts + self.assertEqual(model.default_var_pts, {x: 20}) + + # test geometry + self.assertEqual( + model.default_geometry, {"ECMD particle": {x: {"min": 0, "max": 1}}} + ) + + # test spatial methods + self.assertIsInstance( + model.default_spatial_methods["ECMD particle"], pybamm.FiniteVolume + ) + + # test submesh types + self.assertEqual( + model.default_submesh_types, {"ECMD particle": pybamm.Uniform1DSubMesh} + ) + def test_changing_number_of_rcs(self): options = {"number of rc elements": 0} model = pybamm.equivalent_circuit.Thevenin(options=options) @@ -33,6 +55,11 @@ def test_changing_number_of_rcs(self): model = pybamm.equivalent_circuit.Thevenin(options=options) model.check_well_posedness() + def test_diffusion_element(self): + options = {"diffusion element": "true"} + model = pybamm.equivalent_circuit.Thevenin(options=options) + model.check_well_posedness(post_discretisation=True) + def test_calculate_discharge_energy(self): options = {"calculate discharge energy": "true"} model = pybamm.equivalent_circuit.Thevenin(options=options)