From 0e70d3e94819c7109aa78c4d87e0029a369cc7b9 Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Tue, 27 Aug 2024 10:51:58 +0100 Subject: [PATCH 01/20] Started modifying porosity change submodel to work with phases --- .../porosity/reaction_driven_porosity.py | 71 ++++++++++++------- 1 file changed, 47 insertions(+), 24 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index fc69d0f1fd..3dceb14eba 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -23,38 +23,61 @@ def __init__(self, param, options, x_average): self.x_average = x_average def get_coupled_variables(self, variables): + param = self.param eps_dict = {} for domain in self.options.whole_cell_domains: + delta_eps_k = 0 if domain == "separator": - delta_eps_k = 0 # separator porosity does not change + pass # separator porosity does not change else: - Domain = domain.split()[0].capitalize() - L_sei_k = variables[f"{Domain} total SEI thickness [m]"] - if Domain == "Negative": - L_sei_0 = self.param.n.prim.L_inner_0 + self.param.n.prim.L_outer_0 - elif Domain == "Positive": - L_sei_0 = self.param.p.prim.L_inner_0 + self.param.p.prim.L_outer_0 - L_pl_k = variables[f"{Domain} lithium plating thickness [m]"] - L_dead_k = variables[f"{Domain} dead lithium thickness [m]"] - L_sei_cr_k = variables[f"{Domain} total SEI on cracks thickness [m]"] + dom = domain.split()[0] + Domain = dom.capitalize() roughness_k = variables[f"{Domain} electrode roughness ratio"] + phases = self.options.phases[dom] + for phase in phases: + L_sei_k = variables[f"{Domain} total {phase} SEI thickness [m]"] + if Domain == "Negative": + if phase == "secondary": + L_sei_0 = ( + param.n.sec.L_inner_0 + param.n.sec.L_outer_0 + ) + else: + L_sei_0 = ( + param.n.prim.L_inner_0 + param.n.prim.L_outer_0 + ) + elif Domain == "Positive": + if phase == "secondary": + L_sei_0 = ( + param.p.sec.L_inner_0 + param.p.sec.L_outer_0 + ) + else: + L_sei_0 = ( + param.p.prim.L_inner_0 + param.p.prim.L_outer_0 + ) + L_pl_k = variables[ + f"{Domain} {phase} lithium plating thickness [m]" + ] + L_dead_k = variables[f"{Domain} {phase} dead lithium thickness [m]"] + L_sei_cr_k = variables[ + f"{Domain} total {phase} SEI on cracks thickness [m]" + ] - L_tot = ( - (L_sei_k - L_sei_0) - + L_pl_k - + L_dead_k - + L_sei_cr_k * (roughness_k - 1) - ) + L_tot = ( + (L_sei_k - L_sei_0) + + L_pl_k + + L_dead_k + + L_sei_cr_k * (roughness_k - 1) + ) - a_k = variables[ - f"{Domain} electrode surface area to volume ratio [m-1]" - ] + a_k = variables[ + f"{Domain} electrode {phase} surface area to volume ratio [m-1]" + ] - # This assumes a thin film so curvature effects are neglected. - # They could be included (e.g. for a sphere it is - # a_n * (L_tot + L_tot ** 2 / R_n + L_tot ** # 3 / (3 * R_n ** 2))) - # but it is not clear if it is relevant or not. - delta_eps_k = -a_k * L_tot + # This assumes a thin film so curvature effects are neglected. + # They could be included (e.g. for a sphere it is + # a_n * (L_tot + L_tot ** 2 / R_n + L_tot ** # 3 / (3 * R_n ** 2))) + # but it is not clear if it is relevant or not. + delta_eps_k += -a_k * L_tot domain_param = self.param.domain_params[domain.split()[0]] eps_k = domain_param.epsilon_init + delta_eps_k From b8dd0138d134a9f372276b257ae6f76b626fafa3 Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Tue, 27 Aug 2024 11:03:22 +0100 Subject: [PATCH 02/20] Porosity now works on composite electrode --- .../models/submodels/porosity/reaction_driven_porosity.py | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 3dceb14eba..8a9d195878 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -33,7 +33,11 @@ def get_coupled_variables(self, variables): dom = domain.split()[0] Domain = dom.capitalize() roughness_k = variables[f"{Domain} electrode roughness ratio"] - phases = self.options.phases[dom] + phases_option = getattr(self.options, dom)["particle phases"] + if phases_option == "1": + phases = "" + else: + phases = self.options.phases[dom] for phase in phases: L_sei_k = variables[f"{Domain} total {phase} SEI thickness [m]"] if Domain == "Negative": From e8fc8f47764693a7f98dcee898e7bfcd47ccdbd5 Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Tue, 27 Aug 2024 11:17:23 +0100 Subject: [PATCH 03/20] Now works for non-composite electrode as well! --- .../porosity/reaction_driven_porosity.py | 21 +++++++++++-------- 1 file changed, 12 insertions(+), 9 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 8a9d195878..700d6dc0c4 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -34,12 +34,13 @@ def get_coupled_variables(self, variables): Domain = dom.capitalize() roughness_k = variables[f"{Domain} electrode roughness ratio"] phases_option = getattr(self.options, dom)["particle phases"] - if phases_option == "1": - phases = "" - else: - phases = self.options.phases[dom] + phases = self.options.phases[dom] for phase in phases: - L_sei_k = variables[f"{Domain} total {phase} SEI thickness [m]"] + if phases_option == "1" and phase == "primary": + phase_name = "" + else: + phase_name = phase + " " + L_sei_k = variables[f"{Domain} total {phase_name}SEI thickness [m]"] if Domain == "Negative": if phase == "secondary": L_sei_0 = ( @@ -59,11 +60,13 @@ def get_coupled_variables(self, variables): param.p.prim.L_inner_0 + param.p.prim.L_outer_0 ) L_pl_k = variables[ - f"{Domain} {phase} lithium plating thickness [m]" + f"{Domain} {phase_name}lithium plating thickness [m]" + ] + L_dead_k = variables[ + f"{Domain} {phase_name}dead lithium thickness [m]" ] - L_dead_k = variables[f"{Domain} {phase} dead lithium thickness [m]"] L_sei_cr_k = variables[ - f"{Domain} total {phase} SEI on cracks thickness [m]" + f"{Domain} total {phase_name}SEI on cracks thickness [m]" ] L_tot = ( @@ -74,7 +77,7 @@ def get_coupled_variables(self, variables): ) a_k = variables[ - f"{Domain} electrode {phase} surface area to volume ratio [m-1]" + f"{Domain} electrode {phase_name}surface area to volume ratio [m-1]" ] # This assumes a thin film so curvature effects are neglected. From b76a2ca33d2954f4e3990965b857ce59950c68d9 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 5 Sep 2024 11:04:04 +0000 Subject: [PATCH 04/20] style: pre-commit fixes --- .../porosity/reaction_driven_porosity.py | 16 ++++------------ 1 file changed, 4 insertions(+), 12 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 700d6dc0c4..b9bb898bcb 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -43,22 +43,14 @@ def get_coupled_variables(self, variables): L_sei_k = variables[f"{Domain} total {phase_name}SEI thickness [m]"] if Domain == "Negative": if phase == "secondary": - L_sei_0 = ( - param.n.sec.L_inner_0 + param.n.sec.L_outer_0 - ) + L_sei_0 = param.n.sec.L_inner_0 + param.n.sec.L_outer_0 else: - L_sei_0 = ( - param.n.prim.L_inner_0 + param.n.prim.L_outer_0 - ) + L_sei_0 = param.n.prim.L_inner_0 + param.n.prim.L_outer_0 elif Domain == "Positive": if phase == "secondary": - L_sei_0 = ( - param.p.sec.L_inner_0 + param.p.sec.L_outer_0 - ) + L_sei_0 = param.p.sec.L_inner_0 + param.p.sec.L_outer_0 else: - L_sei_0 = ( - param.p.prim.L_inner_0 + param.p.prim.L_outer_0 - ) + L_sei_0 = param.p.prim.L_inner_0 + param.p.prim.L_outer_0 L_pl_k = variables[ f"{Domain} {phase_name}lithium plating thickness [m]" ] From 6b870412c2cf9a91488b0fee4af1e92a28f70c48 Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Thu, 5 Sep 2024 12:05:48 +0100 Subject: [PATCH 05/20] changelog --- CHANGELOG.md | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/CHANGELOG.md b/CHANGELOG.md index 5203229bd3..4491e93b4c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,5 +1,9 @@ # [Unreleased](https://github.com/pybamm-team/PyBaMM/) +## Features + +- Porosity change now works for composite electrode ([#4417](https://github.com/pybamm-team/PyBaMM/pull/4417)) + ## Optimizations - Removed the `start_step_offset` setting and disabled minimum `dt` warnings for drive cycles with the (`IDAKLUSolver`). ([#4416](https://github.com/pybamm-team/PyBaMM/pull/4416)) From 4758ab88f6ef17c53e094db12a591cba2ae00c48 Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Thu, 5 Sep 2024 12:08:44 +0100 Subject: [PATCH 06/20] style fix --- .../models/submodels/porosity/reaction_driven_porosity.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index b9bb898bcb..2998330c45 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -69,7 +69,8 @@ def get_coupled_variables(self, variables): ) a_k = variables[ - f"{Domain} electrode {phase_name}surface area to volume ratio [m-1]" + f"{Domain} electrode {phase_name}" + "surface area to volume ratio [m-1]" ] # This assumes a thin film so curvature effects are neglected. From adcb65284fe57597e81d8e623987340489568bec Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Thu, 5 Sep 2024 12:52:53 +0100 Subject: [PATCH 07/20] Updated tests --- .../test_lithium_ion/base_lithium_ion_tests.py | 1 + .../test_lithium_ion/base_lithium_ion_tests.py | 2 ++ 2 files changed, 3 insertions(+) diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py index 60e8dfb819..7e7bed7ea5 100644 --- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py +++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py @@ -333,6 +333,7 @@ def test_composite_graphite_silicon_sei(self): "particle phases": ("2", "1"), "open-circuit potential": (("single", "current sigmoid"), "single"), "SEI": "ec reaction limited", + "SEI porosity change": "true", } parameter_values = pybamm.ParameterValues("Chen2020_composite") name = "Negative electrode active material volume fraction" diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py index 9c093c0c65..2772fc56a1 100644 --- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py +++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py @@ -569,6 +569,7 @@ def test_well_posed_composite_different_degradation(self): options = { "particle phases": ("2", "1"), "SEI": ("ec reaction limited", "none"), + "SEI porosity change": "true", "lithium plating": ("reversible", "none"), "open-circuit potential": (("current sigmoid", "single"), "single"), } @@ -577,6 +578,7 @@ def test_well_posed_composite_different_degradation(self): options = { "particle phases": ("2", "1"), "SEI": (("ec reaction limited", "solvent-diffusion limited"), "none"), + "SEI porosity change": "true", "lithium plating": (("reversible", "irreversible"), "none"), "open-circuit potential": (("current sigmoid", "single"), "single"), } From 1985fb169700a88d148f68735fccc69e2f8fd681 Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Sat, 7 Sep 2024 22:52:39 +0100 Subject: [PATCH 08/20] Changed how pref is handled --- .../porosity/reaction_driven_porosity.py | 25 ++++++++++--------- 1 file changed, 13 insertions(+), 12 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 2998330c45..17e93879f0 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -23,7 +23,6 @@ def __init__(self, param, options, x_average): self.x_average = x_average def get_coupled_variables(self, variables): - param = self.param eps_dict = {} for domain in self.options.whole_cell_domains: delta_eps_k = 0 @@ -36,21 +35,23 @@ def get_coupled_variables(self, variables): phases_option = getattr(self.options, dom)["particle phases"] phases = self.options.phases[dom] for phase in phases: - if phases_option == "1" and phase == "primary": + if self.options["particle phases"] == "1": phase_name = "" + pref = "" + elif phases_option == "1" and phase == "primary": + phase_name = "" + pref = "Primary: " else: phase_name = phase + " " + pref = phase.capitalize() + ": " L_sei_k = variables[f"{Domain} total {phase_name}SEI thickness [m]"] - if Domain == "Negative": - if phase == "secondary": - L_sei_0 = param.n.sec.L_inner_0 + param.n.sec.L_outer_0 - else: - L_sei_0 = param.n.prim.L_inner_0 + param.n.prim.L_outer_0 - elif Domain == "Positive": - if phase == "secondary": - L_sei_0 = param.p.sec.L_inner_0 + param.p.sec.L_outer_0 - else: - L_sei_0 = param.p.prim.L_inner_0 + param.p.prim.L_outer_0 + L_inner_0 = pybamm.Parameter( + f"{pref}Initial inner SEI thickness [m]" + ) + L_outer_0 = pybamm.Parameter( + f"{pref}Initial outer SEI thickness [m]" + ) + L_sei_0 = L_inner_0 + L_outer_0 L_pl_k = variables[ f"{Domain} {phase_name}lithium plating thickness [m]" ] From 709d634a62681c11ed4b9be6fe61483f63e2a080 Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Sat, 7 Sep 2024 23:12:36 +0100 Subject: [PATCH 09/20] Added comments to the new if statements --- .../models/submodels/porosity/reaction_driven_porosity.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index c7328936c2..00c4251a72 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -36,12 +36,15 @@ def get_coupled_variables(self, variables): phases = self.options.phases[dom] for phase in phases: if self.options["particle phases"] == "1": + # both electrodes have one phase phase_name = "" pref = "" elif phases_option == "1" and phase == "primary": + # `domain` has one phase phase_name = "" pref = "Primary: " else: + # `domain` has more than one phase phase_name = phase + " " pref = phase.capitalize() + ": " L_sei_k = variables[f"{Domain} total {phase_name}SEI thickness [m]"] From f0a5a7935227beeb6cc466f09cc55a391ff6901c Mon Sep 17 00:00:00 2001 From: Simon O'Kane Date: Mon, 9 Sep 2024 17:28:00 +0100 Subject: [PATCH 10/20] Why did I not think of this before?! --- .../porosity/reaction_driven_porosity.py | 26 +++++++++---------- 1 file changed, 13 insertions(+), 13 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 00c4251a72..d2a01b750e 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -32,29 +32,29 @@ def get_coupled_variables(self, variables): dom = domain.split()[0] Domain = dom.capitalize() roughness_k = variables[f"{Domain} electrode roughness ratio"] + SEI_option = getattr(self.options, dom)["SEI"] phases_option = getattr(self.options, dom)["particle phases"] phases = self.options.phases[dom] for phase in phases: - if self.options["particle phases"] == "1": - # both electrodes have one phase - phase_name = "" - pref = "" - elif phases_option == "1" and phase == "primary": + if phases_option == "1" and phase == "primary": # `domain` has one phase phase_name = "" - pref = "Primary: " + pref = "" else: # `domain` has more than one phase phase_name = phase + " " pref = phase.capitalize() + ": " L_sei_k = variables[f"{Domain} total {phase_name}SEI thickness [m]"] - L_inner_0 = pybamm.Parameter( - f"{pref}Initial inner SEI thickness [m]" - ) - L_outer_0 = pybamm.Parameter( - f"{pref}Initial outer SEI thickness [m]" - ) - L_sei_0 = L_inner_0 + L_outer_0 + if SEI_option == "none": + L_sei_0 = pybamm.Scalar(0) + else: + L_inner_0 = pybamm.Parameter( + f"{pref}Initial inner SEI thickness [m]" + ) + L_outer_0 = pybamm.Parameter( + f"{pref}Initial outer SEI thickness [m]" + ) + L_sei_0 = L_inner_0 + L_outer_0 L_pl_k = variables[ f"{Domain} {phase_name}lithium plating thickness [m]" ] From 029cdeeaf09d49e0002aad24e7f86635c8e22510 Mon Sep 17 00:00:00 2001 From: mohammedasher Date: Fri, 11 Oct 2024 15:35:11 +0100 Subject: [PATCH 11/20] mechanics for composite electrode --- .gitignore | 2 + .../lithium_ion/Chen2020_composite.py | 357 ++++++++++++++++-- .../full_battery_models/base_battery_model.py | 10 +- .../submodels/interface/base_interface.py | 4 +- .../submodels/interface/sei/base_sei.py | 2 +- .../submodels/particle/base_particle.py | 2 +- .../particle_mechanics/base_mechanics.py | 20 +- .../particle_mechanics/crack_propagation.py | 42 +-- .../particle_mechanics/no_mechanics.py | 6 +- .../particle_mechanics/swelling_only.py | 6 +- .../porosity/reaction_driven_porosity.py | 2 +- .../parameters/lithium_ion_parameters.py | 83 ++-- 12 files changed, 404 insertions(+), 132 deletions(-) diff --git a/.gitignore b/.gitignore index 42c76b7c55..7e558a9eb1 100644 --- a/.gitignore +++ b/.gitignore @@ -137,3 +137,5 @@ results/ # tests test_callback.log +PublicPyBaMM_Composite/ +PyBaMM/ diff --git a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py index 69b622a7c5..ece594c788 100644 --- a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py +++ b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py @@ -39,6 +39,81 @@ def graphite_LGM50_electrolyte_exchange_current_density_Chen2020( return m_ref * arrhenius * c_e**0.5 * c_s_surf**0.5 * (c_s_max - c_s_surf) ** 0.5 +def graphite_volume_change_Ai2020(sto, c_s_max): + """ + Graphite particle volume change as a function of stochiometry [1, 2]. + References + ---------- + .. [1] Ai, W., Kraft, L., Sturm, J., Jossen, A., & Wu, B. (2020). + Electrochemical Thermal-Mechanical Modelling of Stress Inhomogeneity in + Lithium-Ion Pouch Cells. Journal of The Electrochemical Society, 167(1), 013512 + DOI: 10.1149/2.0122001JES. + .. [2] Rieger, B., Erhard, S. V., Rumpf, K., & Jossen, A. (2016). + A new method to model the thickness change of a commercial pouch cell + during discharge. Journal of The Electrochemical Society, 163(8), A1566-A1575. + Parameters + ---------- + sto: :class:`pybamm.Symbol` + Electrode stochiometry, dimensionless + should be R-averaged particle concentration + Returns + ------- + t_change:class:`pybamm.Symbol` + volume change, dimensionless, normalised by particle volume + """ + p1 = 145.907 + p2 = -681.229 + p3 = 1334.442 + p4 = -1415.710 + p5 = 873.906 + p6 = -312.528 + p7 = 60.641 + p8 = -5.706 + p9 = 0.386 + p10 = -4.966e-05 + t_change = 50*( + p1 * sto**9 + + p2 * sto**8 + + p3 * sto**7 + + p4 * sto**6 + + p5 * sto**5 + + p6 * sto**4 + + p7 * sto**3 + + p8 * sto**2 + + p9 * sto + + p10 + ) + """ + omega = pybamm.Parameter("Primary: Negative electrode partial molar volume [m3.mol-1]") + t_change = omega * c_s_max * sto + """ + return t_change +def graphite_cracking_rate_Ai2020(T_dim): + """ + Graphite particle cracking rate as a function of temperature [1, 2]. + References + ---------- + .. [1] Ai, W., Kraft, L., Sturm, J., Jossen, A., & Wu, B. (2020). + Electrochemical Thermal-Mechanical Modelling of Stress Inhomogeneity in + Lithium-Ion Pouch Cells. Journal of The Electrochemical Society, 167(1), 013512 + DOI: 10.1149/2.0122001JES. + .. [2] Deshpande, R., Verbrugge, M., Cheng, Y. T., Wang, J., & Liu, P. (2012). + Battery cycle life prediction with coupled chemical degradation and fatigue + mechanics. Journal of the Electrochemical Society, 159(10), A1730. + Parameters + ---------- + T_dim: :class:`pybamm.Symbol` + temperature, [K] + Returns + ------- + k_cr: :class:`pybamm.Symbol` + cracking rate, [m/(Pa.m0.5)^m_cr] + where m_cr is another Paris' law constant + """ + k_cr = 3.9e-21 + Eac_cr = 0 # to be implemented + arrhenius = pybamm.exp(Eac_cr / pybamm.constants.R * (1 / T_dim - 1 / 298.15)) + return k_cr * arrhenius def silicon_ocp_lithiation_Mark2016(sto): """ @@ -128,6 +203,7 @@ def silicon_ocp_delithiation_Mark2016(sto): return U_delithiation + def silicon_LGM50_electrolyte_exchange_current_density_Chen2020( c_e, c_s_surf, c_s_max, T ): @@ -167,6 +243,109 @@ def silicon_LGM50_electrolyte_exchange_current_density_Chen2020( return m_ref * arrhenius * c_e**0.5 * c_s_surf**0.5 * (c_s_max - c_s_surf) ** 0.5 +def silicon_volume_change_Ai2020(sto, c_s_max): + """ + Silicon particle volume change as a function of stochiometry [1, 2]. + References + ---------- + .. [1] Ai, W., Kraft, L., Sturm, J., Jossen, A., & Wu, B. (2020). + Electrochemical Thermal-Mechanical Modelling of Stress Inhomogeneity in + Lithium-Ion Pouch Cells. Journal of The Electrochemical Society, 167(1), 013512 + DOI: 10.1149/2.0122001JES. + .. [2] Rieger, B., Erhard, S. V., Rumpf, K., & Jossen, A. (2016). + A new method to model the thickness change of a commercial pouch cell + during discharge. Journal of The Electrochemical Society, 163(8), A1566-A1575. + Parameters + ---------- + sto: :class:`pybamm.Symbol` + Electrode stochiometry, dimensionless + should be R-averaged particle concentration + Returns + ------- + t_change:class:`pybamm.Symbol` + volume change, dimensionless, normalised by particle volume + """ + p1 = 145.907 + p2 = -681.229 + p3 = 1334.442 + p4 = -1415.710 + p5 = 873.906 + p6 = -312.528 + p7 = 60.641 + p8 = -5.706 + p9 = 0.386 + p10 = -4.966e-05 + t_change = 50*( + p1 * sto**9 + + p2 * sto**8 + + p3 * sto**7 + + p4 * sto**6 + + p5 * sto**5 + + p6 * sto**4 + + p7 * sto**3 + + p8 * sto**2 + + p9 * sto + + p10 + ) + """ + omega = pybamm.Parameter("Secondary: Negative electrode partial molar volume [m3.mol-1]") + t_change = omega * c_s_max * sto + """ + return t_change +def silicon_cracking_rate_Ai2020(T_dim): + """ + Silicon particle cracking rate as a function of temperature [1, 2]. + References + ---------- + .. [1] Ai, W., Kraft, L., Sturm, J., Jossen, A., & Wu, B. (2020). + Electrochemical Thermal-Mechanical Modelling of Stress Inhomogeneity in + Lithium-Ion Pouch Cells. Journal of The Electrochemical Society, 167(1), 013512 + DOI: 10.1149/2.0122001JES. + .. [2] Deshpande, R., Verbrugge, M., Cheng, Y. T., Wang, J., & Liu, P. (2012). + Battery cycle life prediction with coupled chemical degradation and fatigue + mechanics. Journal of the Electrochemical Society, 159(10), A1730. + Parameters + ---------- + T_dim: :class:`pybamm.Symbol` + temperature, [K] + Returns + ------- + k_cr: :class:`pybamm.Symbol` + cracking rate, [m/(Pa.m0.5)^m_cr] + where m_cr is another Paris' law constant + """ + k_cr = 1.0e-22 + Eac_cr = 0 # to be implemented + arrhenius = pybamm.exp(Eac_cr / pybamm.constants.R * (1 / T_dim - 1 / 298.15)) + return k_cr * arrhenius +def nmc_LGM50_ocp_Chen2020(sto): + """ + LG M50 NMC open-circuit potential as a function of stochiometry, fit taken + from [1]. + References + ---------- + .. [1] Chang-Hui Chen, Ferran Brosa Planella, Kieran O’Regan, Dominika Gastol, W. + Dhammika Widanage, and Emma Kendrick. "Development of Experimental Techniques for + Parameterization of Multi-scale Lithium-ion Battery Models." Journal of the + Electrochemical Society 167 (2020): 080534. + Parameters + ---------- + sto: :class:`pybamm.Symbol` + Electrode stochiometry + Returns + ------- + :class:`pybamm.Symbol` + Open-circuit potential + """ + u_eq = ( + -0.8090 * sto + + 4.4875 + - 0.0428 * np.tanh(18.5138 * (sto - 0.5542)) + - 17.7326 * np.tanh(15.7890 * (sto - 0.3117)) + + 17.5842 * np.tanh(15.9308 * (sto - 0.3120)) + ) + return u_eq + def nmc_LGM50_ocp_Chen2020(sto): """ @@ -300,6 +479,58 @@ def electrolyte_conductivity_Nyman2008(c_e, T): return sigma_e +def volume_change_Ai2020(sto, c_s_max): + """ + Particle volume change as a function of stochiometry [1, 2]. + References + ---------- + .. [1] > Ai, W., Kraft, L., Sturm, J., Jossen, A., & Wu, B. (2020). + Electrochemical Thermal-Mechanical Modelling of Stress Inhomogeneity in + Lithium-Ion Pouch Cells. Journal of The Electrochemical Society, 167(1), 013512 + DOI: 10.1149/2.0122001JES. + .. [2] > Rieger, B., Erhard, S. V., Rumpf, K., & Jossen, A. (2016). + A new method to model the thickness change of a commercial pouch cell + during discharge. Journal of The Electrochemical Society, 163(8), A1566-A1575. + Parameters + ---------- + sto: :class:`pybamm.Symbol` + Electrode stochiometry, dimensionless + should be R-averaged particle concentration + Returns + ------- + t_change:class:`pybamm.Symbol` + volume change, dimensionless, normalised by particle volume + """ + omega = pybamm.Parameter("Positive electrode partial molar volume [m3.mol-1]") + t_change = omega * c_s_max * sto + return t_change +def cracking_rate_Ai2020(T_dim): + """ + Particle cracking rate as a function of temperature [1, 2]. + References + ---------- + .. [1] > Ai, W., Kraft, L., Sturm, J., Jossen, A., & Wu, B. (2020). + Electrochemical Thermal-Mechanical Modelling of Stress Inhomogeneity in + Lithium-Ion Pouch Cells. Journal of The Electrochemical Society, 167(1), 013512 + DOI: 10.1149/2.0122001JES. + .. [2] > Deshpande, R., Verbrugge, M., Cheng, Y. T., Wang, J., & Liu, P. (2012). + Battery cycle life prediction with coupled chemical degradation and fatigue + mechanics. Journal of the Electrochemical Society, 159(10), A1730. + Parameters + ---------- + T: :class:`pybamm.Symbol` + temperature, [K] + Returns + ------- + k_cr: :class:`pybamm.Symbol` + cracking rate, [m/(Pa.m0.5)^m_cr] + where m_cr is another Paris' law constant + """ + k_cr = 1.0e-22 + Eac_cr = 0 # to be implemented + arrhenius = pybamm.exp(Eac_cr / pybamm.constants.R * (1 / T_dim - 1 / 298.15)) + return k_cr * arrhenius + # Load data in the appropriate format path, _ = os.path.split(os.path.abspath(__file__)) @@ -327,46 +558,47 @@ def get_parameter_values(): return { "chemistry": "lithium_ion", # sei - "Primary: Ratio of lithium moles to SEI moles": 2.0, - "Primary: Inner SEI reaction proportion": 0.5, - "Primary: Inner SEI partial molar volume [m3.mol-1]": 9.585e-05, - "Primary: Outer SEI partial molar volume [m3.mol-1]": 9.585e-05, - "Primary: SEI reaction exchange current density [A.m-2]": 1.5e-07, - "Primary: SEI resistivity [Ohm.m]": 200000.0, - "Primary: Outer SEI solvent diffusivity [m2.s-1]": 2.5000000000000002e-22, - "Primary: Bulk solvent concentration [mol.m-3]": 2636.0, - "Primary: Inner SEI open-circuit potential [V]": 0.1, - "Primary: Outer SEI open-circuit potential [V]": 0.8, - "Primary: Inner SEI electron conductivity [S.m-1]": 8.95e-14, - "Primary: Inner SEI lithium interstitial diffusivity [m2.s-1]": 1e-20, - "Primary: Lithium interstitial reference concentration [mol.m-3]": 15.0, - "Primary: Initial inner SEI thickness [m]": 2.5e-09, - "Primary: Initial outer SEI thickness [m]": 2.5e-09, - "Primary: EC initial concentration in electrolyte [mol.m-3]": 4541.0, - "Primary: EC diffusivity [m2.s-1]": 2e-18, - "Primary: SEI kinetic rate constant [m.s-1]": 1e-12, - "Primary: SEI open-circuit potential [V]": 0.4, - "Primary: SEI growth activation energy [J.mol-1]": 0.0, - "Secondary: Ratio of lithium moles to SEI moles": 2.0, - "Secondary: Inner SEI reaction proportion": 0.5, - "Secondary: Inner SEI partial molar volume [m3.mol-1]": 9.585e-05, - "Secondary: Outer SEI partial molar volume [m3.mol-1]": 9.585e-05, - "Secondary: SEI reaction exchange current density [A.m-2]": 1.5e-07, - "Secondary: SEI resistivity [Ohm.m]": 200000.0, - "Secondary: Outer SEI solvent diffusivity [m2.s-1]": 2.5000000000000002e-22, - "Secondary: Bulk solvent concentration [mol.m-3]": 2636.0, - "Secondary: Inner SEI open-circuit potential [V]": 0.1, - "Secondary: Outer SEI open-circuit potential [V]": 0.8, - "Secondary: Inner SEI electron conductivity [S.m-1]": 8.95e-14, - "Secondary: Inner SEI lithium interstitial diffusivity [m2.s-1]": 1e-20, - "Secondary: Lithium interstitial reference concentration [mol.m-3]": 15.0, - "Secondary: Initial inner SEI thickness [m]": 2.5e-09, - "Secondary: Initial outer SEI thickness [m]": 2.5e-09, - "Secondary: EC initial concentration in electrolyte [mol.m-3]": 4541.0, - "Secondary: EC diffusivity [m2.s-1]": 2e-18, - "Secondary: SEI kinetic rate constant [m.s-1]": 1e-12, - "Secondary: SEI open-circuit potential [V]": 0.4, - "Secondary: SEI growth activation energy [J.mol-1]": 0.0, + "Primary: Negative ratio of lithium moles to SEI moles": 2.0, + "Primary: Negative inner SEI reaction proportion": 0.5, + "Primary: Negative inner SEI partial molar volume [m3.mol-1]": 9.585e-05, + "Primary: Negative outer SEI partial molar volume [m3.mol-1]": 9.585e-05, + "Primary: Negative SEI reaction exchange current density [A.m-2]": 1.5e-07, + "Primary: Negative SEI resistivity [Ohm.m]": 200000.0, + "Primary: Negative outer SEI solvent diffusivity [m2.s-1]": 2.5000000000000002e-22, + "Primary: Negative bulk solvent concentration [mol.m-3]": 2636.0, + "Primary: Negative inner SEI open-circuit potential [V]": 0.1, + "Primary: Negative outer SEI open-circuit potential [V]": 0.8, + "Primary: Negative inner SEI electron conductivity [S.m-1]": 8.95e-14, + "Primary: Negative inner SEI lithium interstitial diffusivity [m2.s-1]": 1e-20, + "Primary: Negative lithium interstitial reference concentration [mol.m-3]": 15.0, + "Primary: Negative initial inner SEI thickness [m]": 2.5e-09, + "Primary: Negative initial outer SEI thickness [m]": 2.5e-09, + "Primary: Negative EC initial concentration in electrolyte [mol.m-3]": 4541.0, + "Primary: Negative EC diffusivity [m2.s-1]": 2e-18, + "Primary: Negative SEI kinetic rate constant [m.s-1]": 1e-12, + "Primary: Negative SEI open-circuit potential [V]": 0.4, + "Primary: Negative SEI growth activation energy [J.mol-1]": 0.0, + "Secondary: Negative ratio of lithium moles to SEI moles": 2.0, + "Secondary: Negative inner SEI reaction proportion": 0.5, + "Secondary: Negative inner SEI partial molar volume [m3.mol-1]": 9.585e-05, + "Secondary: Negative outer SEI partial molar volume [m3.mol-1]": 9.585e-05, + "Secondary: Negative SEI reaction exchange current density [A.m-2]": 1.5e-07, + "Secondary: Negative SEI resistivity [Ohm.m]": 200000.0, + "Secondary: Negative outer SEI solvent diffusivity [m2.s-1]": 2.5000000000000002e-22, + "Secondary: Negative bulk solvent concentration [mol.m-3]": 2636.0, + "Secondary: Negative inner SEI open-circuit potential [V]": 0.1, + "Secondary: Negative outer SEI open-circuit potential [V]": 0.8, + "Secondary: Negative inner SEI electron conductivity [S.m-1]": 8.95e-14, + "Secondary: Negative inner SEI lithium interstitial diffusivity [m2.s-1]": 1e-20, + "Secondary: Negative lithium interstitial reference concentration [mol.m-3]": 15.0, + "Secondary: Negative initial inner SEI thickness [m]": 2.5e-09, + "Secondary: Negative initial outer SEI thickness [m]": 2.5e-09, + "Secondary: Negative EC initial concentration in electrolyte [mol.m-3]": 4541.0, + "Secondary: Negative EC diffusivity [m2.s-1]": 2e-18, + "Secondary: Negative SEI kinetic rate constant [m.s-1]": 1e-12, + "Secondary: Negative SEI open-circuit potential [V]": 0.4, + "Secondary: Negative SEI growth activation energy [J.mol-1]": 0.0, + "Negative electrode reaction-driven LAM factor [m3.mol-1]": 0.0, "Positive electrode reaction-driven LAM factor [m3.mol-1]": 0.0, # cell "Negative current collector thickness [m]": 1.2e-05, @@ -409,6 +641,22 @@ def get_parameter_values(): "Negative electrode specific heat capacity [J.kg-1.K-1]": 700.0, "Negative electrode thermal conductivity [W.m-1.K-1]": 1.7, "Primary: Negative electrode OCP entropic change [V.K-1]": 0.0, + "Primary: Negative electrode Poisson's ratio": 0.3, + "Primary: Negative electrode Young's modulus [Pa]": 15000000000.0, + "Primary: Negative electrode reference concentration for free of deformation [mol.m-3]" + "": 0.0, + "Primary: Negative electrode partial molar volume [m3.mol-1]": 3.1e-06, + "Primary: Negative electrode volume change": graphite_volume_change_Ai2020, + "Primary: Negative electrode initial crack length [m]": 2e-08, + "Primary: Negative electrode initial crack width [m]": 1.5e-08, + "Primary: Negative electrode number of cracks per unit area [m-2]": 3180000000000000.0, + "Primary: Negative electrode Paris' law constant b": 1.12, + "Primary: Negative electrode Paris' law constant m": 2.2, + "Primary: Negative electrode cracking rate": graphite_cracking_rate_Ai2020, + "Negative electrode activation energy for cracking rate [J.mol-1]": 0, + "Primary: Negative electrode LAM constant proportional term [s-1]": 2.7778e-07, + "Primary: Negative electrode LAM constant exponential term": 2.0, + "Primary: Negative electrode critical stress [Pa]": 60000000.0, "Secondary: Maximum concentration in negative electrode [mol.m-3]": 278000.0, "Secondary: Initial concentration in negative electrode [mol.m-3]": 276610.0, "Secondary: Negative particle diffusivity [m2.s-1]": 1.67e-14, @@ -422,6 +670,20 @@ def get_parameter_values(): "": silicon_LGM50_electrolyte_exchange_current_density_Chen2020, "Secondary: Negative electrode density [kg.m-3]": 2650.0, "Secondary: Negative electrode OCP entropic change [V.K-1]": 0.0, + "Secondary: Negative electrode Poisson's ratio": 0.22, + "Secondary: Negative electrode Young's modulus [Pa]": 50000000000.0, + "Secondary: Negative electrode reference concentration for free of deformation [mol.m-3]": 0, + "Secondary: Negative electrode partial molar volume [m3.mol-1]": 1.20e-05, + "Secondary: Negative electrode volume change": silicon_volume_change_Ai2020, + "Secondary: Negative electrode initial crack length [m]": 2e-08, + "Secondary: Negative electrode initial crack width [m]": 1.5e-08, + "Secondary: Negative electrode number of cracks per unit area [m-2]": 3180000000000000.0, + "Secondary: Negative electrode Paris' law constant b": 1.12, + "Secondary: Negative electrode Paris' law constant m": 2.2, + "Secondary: Negative electrode cracking rate": silicon_cracking_rate_Ai2020, + "Secondary: Negative electrode LAM constant proportional term [s-1]": 2.7778e-07, + "Secondary: Negative electrode LAM constant exponential term": 2.0, + "Secondary: Negative electrode critical stress [Pa]": 720000000, # positive electrode "Positive electrode conductivity [S.m-1]": 0.18, "Maximum concentration in positive electrode [mol.m-3]": 63104.0, @@ -440,6 +702,21 @@ def get_parameter_values(): "Positive electrode specific heat capacity [J.kg-1.K-1]": 700.0, "Positive electrode thermal conductivity [W.m-1.K-1]": 2.1, "Positive electrode OCP entropic change [V.K-1]": 0.0, + "Positive electrode Poisson's ratio": 0.2, + "Positive electrode Young's modulus [Pa]": 375000000000.0, + "Positive electrode reference concentration for free of deformation [mol.m-3]" + "": 0.0, + "Positive electrode partial molar volume [m3.mol-1]": 1.25e-05, + "Positive electrode volume change": volume_change_Ai2020, + "Positive electrode initial crack length [m]": 2e-08, + "Positive electrode initial crack width [m]": 1.5e-08, + "Positive electrode number of cracks per unit area [m-2]": 3180000000000000.0, + "Positive electrode Paris' law constant b": 1.12, + "Positive electrode Paris' law constant m": 2.2, + "Positive electrode cracking rate": cracking_rate_Ai2020, + "Positive electrode LAM constant proportional term [s-1]": 2.7778e-07, + "Positive electrode LAM constant exponential term": 2.0, + "Positive electrode critical stress [Pa]": 375000000.0, # separator "Separator porosity": 0.47, "Separator Bruggeman coefficient (electrolyte)": 1.5, diff --git a/src/pybamm/models/full_battery_models/base_battery_model.py b/src/pybamm/models/full_battery_models/base_battery_model.py index 5340d685e3..b9fbc6228f 100644 --- a/src/pybamm/models/full_battery_models/base_battery_model.py +++ b/src/pybamm/models/full_battery_models/base_battery_model.py @@ -594,12 +594,12 @@ def __init__(self, extra_options): if isinstance(options["stress-induced diffusion"], str): if ( options["stress-induced diffusion"] == "true" - and options["particle mechanics"] == "none" + #and options["particle mechanics"] == "none" ): - raise pybamm.OptionError( - "cannot have stress-induced diffusion without a particle " - "mechanics model" - ) + # raise pybamm.OptionError( + # "cannot have stress-induced diffusion without a particle " + # "mechanics model" + # ) if options["working electrode"] != "both": if options["thermal"] == "x-full": diff --git a/src/pybamm/models/submodels/interface/base_interface.py b/src/pybamm/models/submodels/interface/base_interface.py index 0ad08d5454..062ddc2d27 100644 --- a/src/pybamm/models/submodels/interface/base_interface.py +++ b/src/pybamm/models/submodels/interface/base_interface.py @@ -306,9 +306,9 @@ def _get_standard_volumetric_current_density_variables(self, variables): a_j_av = pybamm.x_average(a_j) if reaction_name == "SEI on cracks ": - roughness = variables[f"{Domain} electrode roughness ratio"] - 1 + roughness = variables[f"{Domain} electrode {self.phase_name}roughness ratio"] - 1 roughness_av = ( - variables[f"X-averaged {domain} electrode roughness ratio"] - 1 + variables[f"X-averaged {domain} electrode {self.phase_name}roughness ratio"] - 1 ) else: roughness = 1 diff --git a/src/pybamm/models/submodels/interface/sei/base_sei.py b/src/pybamm/models/submodels/interface/sei/base_sei.py index 6caeac887d..1cf2d06c09 100644 --- a/src/pybamm/models/submodels/interface/sei/base_sei.py +++ b/src/pybamm/models/submodels/interface/sei/base_sei.py @@ -218,7 +218,7 @@ def _get_standard_concentration_variables(self, variables): elif self.reaction == "SEI on cracks": L_inner_cr = variables[f"{Domain} inner {reaction_name}thickness [m]"] L_outer_cr = variables[f"{Domain} outer {reaction_name}thickness [m]"] - roughness = variables[f"{Domain} electrode roughness ratio"] + roughness = variables[f"{Domain} electrode {self.phase_name}roughness ratio"] n_inner_cr = L_inner_cr * L_to_n_inner * (roughness - 1) n_outer_cr = L_outer_cr * L_to_n_outer * (roughness - 1) diff --git a/src/pybamm/models/submodels/particle/base_particle.py b/src/pybamm/models/submodels/particle/base_particle.py index b774e58a0c..6efa707ee3 100644 --- a/src/pybamm/models/submodels/particle/base_particle.py +++ b/src/pybamm/models/submodels/particle/base_particle.py @@ -60,7 +60,7 @@ def _get_effective_diffusivity(self, c, T, current): E = pybamm.r_average(phase_param.E(sto, T)) nu = phase_param.nu theta_M = Omega / (self.param.R * T) * (2 * Omega * E) / (9 * (1 - nu)) - stress_factor = 1 + theta_M * (c - domain_param.c_0) + stress_factor = 1 + theta_M * (c - phase_param.c_0) else: stress_factor = 1 diff --git a/src/pybamm/models/submodels/particle_mechanics/base_mechanics.py b/src/pybamm/models/submodels/particle_mechanics/base_mechanics.py index 1301722da0..02093a753b 100644 --- a/src/pybamm/models/submodels/particle_mechanics/base_mechanics.py +++ b/src/pybamm/models/submodels/particle_mechanics/base_mechanics.py @@ -23,15 +23,15 @@ class BaseMechanics(pybamm.BaseSubModel): """ - def __init__(self, param, domain, options, phase="primary"): + def __init__(self, param, domain, options, phase): super().__init__(param, domain, options=options, phase=phase) def _get_standard_variables(self, l_cr): domain, Domain = self.domain_Domain l_cr_av = pybamm.x_average(l_cr) variables = { - f"{Domain} particle crack length [m]": l_cr, - f"X-averaged {domain} particle crack length [m]": l_cr_av, + f"{Domain} {self.phase_name}particle crack length [m]": l_cr, + f"X-averaged {domain} {self.phase_name}particle crack length [m]": l_cr_av, } return variables @@ -61,7 +61,7 @@ def _get_mechanical_results(self, variables): sto = variables[f"{Domain} {phase_name}particle concentration"] Omega = pybamm.r_average(phase_param.Omega(sto, T)) R0 = phase_param.R - c_0 = domain_param.c_0 + c_0 = phase_param.c_0 E0 = pybamm.r_average(phase_param.E(sto, T)) nu = phase_param.nu L0 = domain_param.L @@ -129,19 +129,19 @@ def _get_standard_surface_variables(self, variables): domain, Domain = self.domain_Domain phase_name = self.phase_name - l_cr = variables[f"{Domain} particle crack length [m]"] + l_cr = variables[f"{Domain} {self.phase_name}particle crack length [m]"] a = variables[ f"{Domain} electrode {phase_name}surface area to volume ratio [m-1]" ] - rho_cr = self.domain_param.rho_cr - w_cr = self.domain_param.w_cr + rho_cr = self.phase_param.rho_cr + w_cr = self.phase_param.w_cr roughness = 1 + 2 * l_cr * rho_cr * w_cr # ratio of cracks to normal surface a_cr = (roughness - 1) * a # crack surface area to volume ratio roughness_xavg = pybamm.x_average(roughness) variables = { - f"{Domain} crack surface to volume ratio [m-1]": a_cr, - f"{Domain} electrode roughness ratio": roughness, - f"X-averaged {domain} electrode roughness ratio": roughness_xavg, + f"{Domain} {self.phase_name}crack surface to volume ratio [m-1]": a_cr, + f"{Domain} electrode {self.phase_name}roughness ratio": roughness, + f"X-averaged {domain} electrode {self.phase_name}roughness ratio": roughness_xavg, } return variables diff --git a/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py b/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py index b51f1d1ebd..3aa393bcc7 100644 --- a/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py +++ b/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py @@ -27,7 +27,7 @@ class CrackPropagation(BaseMechanics): """ - def __init__(self, param, domain, x_average, options, phase="primary"): + def __init__(self, param, domain, x_average, options, phase): super().__init__(param, domain, options, phase) self.x_average = x_average @@ -39,17 +39,17 @@ def get_fundamental_variables(self): if self.x_average is True: l_cr_av = pybamm.Variable( - f"X-averaged {domain} particle crack length [m]", + f"X-averaged {domain} {self.phase_name}particle crack length [m]", domain="current collector", - scale=self.domain_param.l_cr_0, + scale=self.phase_param.l_cr_0, ) l_cr = pybamm.PrimaryBroadcast(l_cr_av, f"{domain} electrode") else: l_cr = pybamm.Variable( - f"{Domain} particle crack length [m]", + f"{Domain} {self.phase_name}particle crack length [m]", domain=f"{domain} electrode", auxiliary_domains={"secondary": "current collector"}, - scale=self.domain_param.l_cr_0, + scale=self.phase_param.l_cr_0, ) variables = self._get_standard_variables(l_cr) @@ -62,18 +62,18 @@ def get_coupled_variables(self, variables): variables.update(self._get_standard_surface_variables(variables)) variables.update(self._get_mechanical_results(variables)) T = variables[f"{Domain} electrode temperature [K]"] - k_cr = self.domain_param.k_cr(T) - m_cr = self.domain_param.m_cr - b_cr = self.domain_param.b_cr - stress_t_surf = variables[f"{Domain} particle surface tangential stress [Pa]"] - l_cr = variables[f"{Domain} particle crack length [m]"] + k_cr = self.phase_param.k_cr(T) + m_cr = self.phase_param.m_cr + b_cr = self.phase_param.b_cr + stress_t_surf = variables[f"{Domain} {self.phase_name}particle surface tangential stress [Pa]"] + l_cr = variables[f"{Domain} {self.phase_name}particle crack length [m]"] # # compressive stress will not lead to crack propagation dK_SIF = stress_t_surf * b_cr * pybamm.sqrt(np.pi * l_cr) * (stress_t_surf >= 0) dl_cr = k_cr * (dK_SIF**m_cr) / 3600 # divide by 3600 to replace t0_cr variables.update( { - f"{Domain} particle cracking rate [m.s-1]": dl_cr, - f"X-averaged {domain} particle cracking rate [m.s-1]": pybamm.x_average( + f"{Domain} {self.phase_name}particle cracking rate [m.s-1]": dl_cr, + f"X-averaged {domain} {self.phase_name}particle cracking rate [m.s-1]": pybamm.x_average( dl_cr ), } @@ -84,21 +84,21 @@ def set_rhs(self, variables): domain, Domain = self.domain_Domain if self.x_average is True: - l_cr = variables[f"X-averaged {domain} particle crack length [m]"] - dl_cr = variables[f"X-averaged {domain} particle cracking rate [m.s-1]"] + l_cr = variables[f"X-averaged {domain} {self.phase_name}particle crack length [m]"] + dl_cr = variables[f"X-averaged {domain} {self.phase_name}particle cracking rate [m.s-1]"] else: - l_cr = variables[f"{Domain} particle crack length [m]"] - dl_cr = variables[f"{Domain} particle cracking rate [m.s-1]"] + l_cr = variables[f"{Domain} {self.phase_name}particle crack length [m]"] + dl_cr = variables[f"{Domain} {self.phase_name}particle cracking rate [m.s-1]"] self.rhs = {l_cr: dl_cr} def set_initial_conditions(self, variables): domain, Domain = self.domain_Domain - l_cr_0 = self.domain_param.l_cr_0 + l_cr_0 = self.phase_param.l_cr_0 if self.x_average is True: l_cr = variables[f"X-averaged {domain} particle crack length [m]"] else: - l_cr = variables[f"{Domain} particle crack length [m]"] + l_cr = variables[f"X-averaged {domain} {self.phase_name}particle crack length [m]"] l_cr_0 = pybamm.PrimaryBroadcast(l_cr_0, f"{domain} electrode") self.initial_conditions = {l_cr: l_cr_0} @@ -108,10 +108,10 @@ def add_events_from(self, variables): if self.x_average is True: l_cr = variables[f"X-averaged {domain} particle crack length [m]"] else: - l_cr = variables[f"{Domain} particle crack length [m]"] + l_cr = variables[f"X-averaged {domain} {self.phase_name}particle crack length [m]"] self.events.append( pybamm.Event( - f"{domain} particle crack length larger than particle radius", - 1 - pybamm.max(l_cr) / self.domain_param.prim.R_typ, + f"{domain} {self.phase_name}particle crack length larger than particle radius", + 1 - pybamm.max(l_cr) / self.phase_param.R_typ, ) ) diff --git a/src/pybamm/models/submodels/particle_mechanics/no_mechanics.py b/src/pybamm/models/submodels/particle_mechanics/no_mechanics.py index 1746d4781f..9359c9e522 100644 --- a/src/pybamm/models/submodels/particle_mechanics/no_mechanics.py +++ b/src/pybamm/models/submodels/particle_mechanics/no_mechanics.py @@ -22,7 +22,7 @@ class NoMechanics(BaseMechanics): Phase of the particle (default is "primary") """ - def __init__(self, param, domain, options, phase="primary"): + def __init__(self, param, domain, options, phase): super().__init__(param, domain, options, phase) def get_fundamental_variables(self): @@ -35,8 +35,8 @@ def get_fundamental_variables(self): variables = self._get_standard_variables(zero) variables.update( { - f"{Domain} particle cracking rate [m.s-1]": zero, - f"X-averaged {domain} particle cracking rate [m.s-1]": zero_av, + f"{Domain} {self.phase_name}particle cracking rate [m.s-1]": zero, + f"X-averaged {domain} {self.phase_name}particle cracking rate [m.s-1]": zero_av, } ) return variables diff --git a/src/pybamm/models/submodels/particle_mechanics/swelling_only.py b/src/pybamm/models/submodels/particle_mechanics/swelling_only.py index 6dd6fde489..d445de43d4 100644 --- a/src/pybamm/models/submodels/particle_mechanics/swelling_only.py +++ b/src/pybamm/models/submodels/particle_mechanics/swelling_only.py @@ -23,7 +23,7 @@ class SwellingOnly(BaseMechanics): """ - def __init__(self, param, domain, options, phase="primary"): + def __init__(self, param, domain, options, phase): super().__init__(param, domain, options, phase) pybamm.citations.register("Ai2019") @@ -39,8 +39,8 @@ def get_fundamental_variables(self): variables = self._get_standard_variables(zero) variables.update( { - f"{Domain} particle cracking rate [m.s-1]": zero, - f"X-averaged {domain} particle cracking rate [m.s-1]": zero_av, + f"{Domain} {self.phase_name}particle cracking rate [m.s-1]": zero, + f"X-averaged {domain} {self.phase_name}particle cracking rate [m.s-1]": zero_av, } ) return variables diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 989c90cd9c..1f297f6b1b 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -37,7 +37,7 @@ def get_coupled_variables(self, variables): L_pl_k = variables[f"{Domain} lithium plating thickness [m]"] L_dead_k = variables[f"{Domain} dead lithium thickness [m]"] L_sei_cr_k = variables[f"{Domain} total SEI on cracks thickness [m]"] - roughness_k = variables[f"{Domain} electrode roughness ratio"] + roughness_k = variables[f"{Domain} electrode {self.phase_name}roughness ratio"] L_tot = ( (L_sei_k - L_sei_0) diff --git a/src/pybamm/parameters/lithium_ion_parameters.py b/src/pybamm/parameters/lithium_ion_parameters.py index 3902242d78..884e5ef1f3 100644 --- a/src/pybamm/parameters/lithium_ion_parameters.py +++ b/src/pybamm/parameters/lithium_ion_parameters.py @@ -78,6 +78,7 @@ def _set_parameters(self): self.voltage_high_cut = self.elec.voltage_high_cut self.ocp_soc_0 = self.elec.ocp_soc_0 self.ocp_soc_100 = self.elec.ocp_soc_100 + # Domain parameters for domain in self.domain_params.values(): @@ -268,20 +269,7 @@ def _set_parameters(self): self.b_s = self.geo.b_s self.tau_s = self.geo.tau_s - # Mechanical parameters - self.c_0 = pybamm.Parameter( - f"{Domain} electrode reference concentration for free of deformation " - "[mol.m-3]" - ) - - self.l_cr_0 = pybamm.Parameter(f"{Domain} electrode initial crack length [m]") - self.w_cr = pybamm.Parameter(f"{Domain} electrode initial crack width [m]") - self.rho_cr = pybamm.Parameter( - f"{Domain} electrode number of cracks per unit area [m-2]" - ) - self.b_cr = pybamm.Parameter(f"{Domain} electrode Paris' law constant b") - self.m_cr = pybamm.Parameter(f"{Domain} electrode Paris' law constant m") - + # Utilisation parameters self.u_init = pybamm.Parameter( f"Initial {domain} electrode interface utilisation" @@ -306,14 +294,6 @@ def sigma(self, T): f"{Domain} electrode conductivity [S.m-1]", inputs ) - def k_cr(self, T): - """ - Cracking rate for the electrode; - """ - Domain = self.domain.capitalize() - return pybamm.FunctionParameter( - f"{Domain} electrode cracking rate", {"Temperature [K]": T} - ) def LAM_rate_current(self, i, T): """ @@ -507,8 +487,19 @@ def _set_parameters(self): if self.options["particle shape"] == "spherical": self.a_typ = 3 * pybamm.xyz_average(self.epsilon_s) / self.R_typ - # Mechanical property + # Mechanical parameters self.nu = pybamm.Parameter(f"{pref}{Domain} electrode Poisson's ratio") + self.c_0 = pybamm.Parameter( + f"{pref}{Domain} electrode reference concentration for free of deformation " + "[mol.m-3]" + ) + self.l_cr_0 = pybamm.Parameter(f"{pref}{Domain} electrode initial crack length [m]") + self.w_cr = pybamm.Parameter(f"{pref}{Domain} electrode initial crack width [m]") + self.rho_cr = pybamm.Parameter( + f"{pref}{Domain} electrode number of cracks per unit area [m-2]" + ) + self.b_cr = pybamm.Parameter(f"{pref}{Domain} electrode Paris' law constant b") + self.m_cr = pybamm.Parameter(f"{pref}{Domain} electrode Paris' law constant m") # Loss of active material parameters self.m_LAM = pybamm.Parameter( @@ -774,6 +765,29 @@ def dxdU(self, U, T): dxdU += self.dxdU_j(U, T, i) return dxdU + def Omega(self, sto, T): + """Dimensional partial molar volume of Li in solid solution [m3.mol-1]""" + Domain = self.domain.capitalize() + inputs = {f"{self.phase_prefactor}{Domain} particle stoichiometry": sto, "Temperature [K]": T} + return pybamm.FunctionParameter( + f"{self.phase_prefactor}{Domain} electrode partial molar volume [m3.mol-1]", inputs + ) + def E(self, sto, T): + """Dimensional Young's modulus""" + Domain = self.domain.capitalize() + inputs = {f"{self.phase_prefactor}{Domain} particle stoichiometry": sto, "Temperature [K]": T} + return pybamm.FunctionParameter( + f"{self.phase_prefactor}{Domain} electrode Young's modulus [Pa]", inputs + ) + def k_cr(self, T): + """ + Cracking rate for the electrode; + """ + Domain = self.domain.capitalize() + return pybamm.FunctionParameter( + f"{self.phase_prefactor}{Domain} electrode cracking rate", {"Temperature [K]": T} + ) + def t_change(self, sto): """ Volume change for the electrode; sto should be R-averaged @@ -782,29 +796,8 @@ def t_change(self, sto): return pybamm.FunctionParameter( f"{self.phase_prefactor}{Domain} electrode volume change", { - f"{Domain} particle stoichiometry": sto, + f"{self.phase_prefactor}{Domain} electrode volume change", }, ) - def Omega(self, sto, T): - """Dimensional partial molar volume of Li in solid solution [m3.mol-1]""" - domain, Domain = self.domain_Domain - inputs = { - f"{self.phase_prefactor} particle stoichiometry": sto, - "Temperature [K]": T, - } - return pybamm.FunctionParameter( - f"{self.phase_prefactor}{Domain} electrode partial molar volume [m3.mol-1]", - inputs, - ) - def E(self, sto, T): - """Dimensional Young's modulus""" - domain, Domain = self.domain_Domain - inputs = { - f"{self.phase_prefactor} particle stoichiometry": sto, - "Temperature [K]": T, - } - return pybamm.FunctionParameter( - f"{self.phase_prefactor}{Domain} electrode Young's modulus [Pa]", inputs - ) From b9a02bbc589bda7f8be0b273cc8a49186b2041e8 Mon Sep 17 00:00:00 2001 From: mohammedasher Date: Tue, 15 Oct 2024 14:25:12 +0100 Subject: [PATCH 12/20] Fix style issues: remove duplicate function, fix indentation, and remove unused variable --- .../lithium_ion/Chen2020_composite.py | 33 ------------------- .../full_battery_models/base_battery_model.py | 1 + .../submodels/particle/base_particle.py | 1 - 3 files changed, 1 insertion(+), 34 deletions(-) diff --git a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py index ece594c788..4406f6044c 100644 --- a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py +++ b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py @@ -347,39 +347,6 @@ def nmc_LGM50_ocp_Chen2020(sto): return u_eq -def nmc_LGM50_ocp_Chen2020(sto): - """ - LG M50 NMC open-circuit potential as a function of stoichiometry, fit taken - from [1]. - - References - ---------- - .. [1] Chang-Hui Chen, Ferran Brosa Planella, Kieran O’Regan, Dominika Gastol, W. - Dhammika Widanage, and Emma Kendrick. "Development of Experimental Techniques for - Parameterization of Multi-scale Lithium-ion Battery Models." Journal of the - Electrochemical Society 167 (2020): 080534. - - Parameters - ---------- - sto: :class:`pybamm.Symbol` - Electrode stoichiometry - - Returns - ------- - :class:`pybamm.Symbol` - Open-circuit potential - """ - - u_eq = ( - -0.8090 * sto - + 4.4875 - - 0.0428 * np.tanh(18.5138 * (sto - 0.5542)) - - 17.7326 * np.tanh(15.7890 * (sto - 0.3117)) - + 17.5842 * np.tanh(15.9308 * (sto - 0.3120)) - ) - - return u_eq - def nmc_LGM50_electrolyte_exchange_current_density_Chen2020(c_e, c_s_surf, c_s_max, T): """ diff --git a/src/pybamm/models/full_battery_models/base_battery_model.py b/src/pybamm/models/full_battery_models/base_battery_model.py index b9fbc6228f..8bf8b2acaf 100644 --- a/src/pybamm/models/full_battery_models/base_battery_model.py +++ b/src/pybamm/models/full_battery_models/base_battery_model.py @@ -596,6 +596,7 @@ def __init__(self, extra_options): options["stress-induced diffusion"] == "true" #and options["particle mechanics"] == "none" ): + pass # raise pybamm.OptionError( # "cannot have stress-induced diffusion without a particle " # "mechanics model" diff --git a/src/pybamm/models/submodels/particle/base_particle.py b/src/pybamm/models/submodels/particle/base_particle.py index 6efa707ee3..267bf9cd30 100644 --- a/src/pybamm/models/submodels/particle/base_particle.py +++ b/src/pybamm/models/submodels/particle/base_particle.py @@ -29,7 +29,6 @@ def __init__(self, param, domain, options, phase="primary"): def _get_effective_diffusivity(self, c, T, current): domain, Domain = self.domain_Domain - domain_param = self.domain_param phase_param = self.phase_param domain_options = getattr(self.options, domain) From 1056b7c3cdc1e0675ff83a54ee99314d69afabe7 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 15 Oct 2024 13:26:06 +0000 Subject: [PATCH 13/20] style: pre-commit fixes --- .../lithium_ion/Chen2020_composite.py | 18 +++++++-- .../full_battery_models/base_battery_model.py | 2 +- .../submodels/interface/base_interface.py | 9 ++++- .../submodels/interface/sei/base_sei.py | 4 +- .../particle_mechanics/crack_propagation.py | 28 ++++++++++---- .../porosity/reaction_driven_porosity.py | 4 +- .../parameters/lithium_ion_parameters.py | 37 ++++++++++++------- 7 files changed, 71 insertions(+), 31 deletions(-) diff --git a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py index 4406f6044c..ff5722f7e2 100644 --- a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py +++ b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py @@ -39,6 +39,7 @@ def graphite_LGM50_electrolyte_exchange_current_density_Chen2020( return m_ref * arrhenius * c_e**0.5 * c_s_surf**0.5 * (c_s_max - c_s_surf) ** 0.5 + def graphite_volume_change_Ai2020(sto, c_s_max): """ Graphite particle volume change as a function of stochiometry [1, 2]. @@ -71,7 +72,7 @@ def graphite_volume_change_Ai2020(sto, c_s_max): p8 = -5.706 p9 = 0.386 p10 = -4.966e-05 - t_change = 50*( + t_change = 50 * ( p1 * sto**9 + p2 * sto**8 + p3 * sto**7 @@ -88,6 +89,8 @@ def graphite_volume_change_Ai2020(sto, c_s_max): t_change = omega * c_s_max * sto """ return t_change + + def graphite_cracking_rate_Ai2020(T_dim): """ Graphite particle cracking rate as a function of temperature [1, 2]. @@ -115,6 +118,7 @@ def graphite_cracking_rate_Ai2020(T_dim): arrhenius = pybamm.exp(Eac_cr / pybamm.constants.R * (1 / T_dim - 1 / 298.15)) return k_cr * arrhenius + def silicon_ocp_lithiation_Mark2016(sto): """ silicon Open-circuit Potential (OCP) as a a function of the @@ -203,7 +207,6 @@ def silicon_ocp_delithiation_Mark2016(sto): return U_delithiation - def silicon_LGM50_electrolyte_exchange_current_density_Chen2020( c_e, c_s_surf, c_s_max, T ): @@ -243,6 +246,7 @@ def silicon_LGM50_electrolyte_exchange_current_density_Chen2020( return m_ref * arrhenius * c_e**0.5 * c_s_surf**0.5 * (c_s_max - c_s_surf) ** 0.5 + def silicon_volume_change_Ai2020(sto, c_s_max): """ Silicon particle volume change as a function of stochiometry [1, 2]. @@ -275,7 +279,7 @@ def silicon_volume_change_Ai2020(sto, c_s_max): p8 = -5.706 p9 = 0.386 p10 = -4.966e-05 - t_change = 50*( + t_change = 50 * ( p1 * sto**9 + p2 * sto**8 + p3 * sto**7 @@ -292,6 +296,8 @@ def silicon_volume_change_Ai2020(sto, c_s_max): t_change = omega * c_s_max * sto """ return t_change + + def silicon_cracking_rate_Ai2020(T_dim): """ Silicon particle cracking rate as a function of temperature [1, 2]. @@ -318,6 +324,8 @@ def silicon_cracking_rate_Ai2020(T_dim): Eac_cr = 0 # to be implemented arrhenius = pybamm.exp(Eac_cr / pybamm.constants.R * (1 / T_dim - 1 / 298.15)) return k_cr * arrhenius + + def nmc_LGM50_ocp_Chen2020(sto): """ LG M50 NMC open-circuit potential as a function of stochiometry, fit taken @@ -347,7 +355,6 @@ def nmc_LGM50_ocp_Chen2020(sto): return u_eq - def nmc_LGM50_electrolyte_exchange_current_density_Chen2020(c_e, c_s_surf, c_s_max, T): """ Exchange-current density for Butler-Volmer reactions between NMC and LiPF6 in @@ -446,6 +453,7 @@ def electrolyte_conductivity_Nyman2008(c_e, T): return sigma_e + def volume_change_Ai2020(sto, c_s_max): """ Particle volume change as a function of stochiometry [1, 2]. @@ -471,6 +479,8 @@ def volume_change_Ai2020(sto, c_s_max): omega = pybamm.Parameter("Positive electrode partial molar volume [m3.mol-1]") t_change = omega * c_s_max * sto return t_change + + def cracking_rate_Ai2020(T_dim): """ Particle cracking rate as a function of temperature [1, 2]. diff --git a/src/pybamm/models/full_battery_models/base_battery_model.py b/src/pybamm/models/full_battery_models/base_battery_model.py index 8bf8b2acaf..63189c0d78 100644 --- a/src/pybamm/models/full_battery_models/base_battery_model.py +++ b/src/pybamm/models/full_battery_models/base_battery_model.py @@ -594,7 +594,7 @@ def __init__(self, extra_options): if isinstance(options["stress-induced diffusion"], str): if ( options["stress-induced diffusion"] == "true" - #and options["particle mechanics"] == "none" + # and options["particle mechanics"] == "none" ): pass # raise pybamm.OptionError( diff --git a/src/pybamm/models/submodels/interface/base_interface.py b/src/pybamm/models/submodels/interface/base_interface.py index 062ddc2d27..b613a82941 100644 --- a/src/pybamm/models/submodels/interface/base_interface.py +++ b/src/pybamm/models/submodels/interface/base_interface.py @@ -306,9 +306,14 @@ def _get_standard_volumetric_current_density_variables(self, variables): a_j_av = pybamm.x_average(a_j) if reaction_name == "SEI on cracks ": - roughness = variables[f"{Domain} electrode {self.phase_name}roughness ratio"] - 1 + roughness = ( + variables[f"{Domain} electrode {self.phase_name}roughness ratio"] - 1 + ) roughness_av = ( - variables[f"X-averaged {domain} electrode {self.phase_name}roughness ratio"] - 1 + variables[ + f"X-averaged {domain} electrode {self.phase_name}roughness ratio" + ] + - 1 ) else: roughness = 1 diff --git a/src/pybamm/models/submodels/interface/sei/base_sei.py b/src/pybamm/models/submodels/interface/sei/base_sei.py index 1cf2d06c09..df17cf3ac7 100644 --- a/src/pybamm/models/submodels/interface/sei/base_sei.py +++ b/src/pybamm/models/submodels/interface/sei/base_sei.py @@ -218,7 +218,9 @@ def _get_standard_concentration_variables(self, variables): elif self.reaction == "SEI on cracks": L_inner_cr = variables[f"{Domain} inner {reaction_name}thickness [m]"] L_outer_cr = variables[f"{Domain} outer {reaction_name}thickness [m]"] - roughness = variables[f"{Domain} electrode {self.phase_name}roughness ratio"] + roughness = variables[ + f"{Domain} electrode {self.phase_name}roughness ratio" + ] n_inner_cr = L_inner_cr * L_to_n_inner * (roughness - 1) n_outer_cr = L_outer_cr * L_to_n_outer * (roughness - 1) diff --git a/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py b/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py index 3aa393bcc7..c16e4e0b96 100644 --- a/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py +++ b/src/pybamm/models/submodels/particle_mechanics/crack_propagation.py @@ -46,7 +46,7 @@ def get_fundamental_variables(self): l_cr = pybamm.PrimaryBroadcast(l_cr_av, f"{domain} electrode") else: l_cr = pybamm.Variable( - f"{Domain} {self.phase_name}particle crack length [m]", + f"{Domain} {self.phase_name}particle crack length [m]", domain=f"{domain} electrode", auxiliary_domains={"secondary": "current collector"}, scale=self.phase_param.l_cr_0, @@ -65,7 +65,9 @@ def get_coupled_variables(self, variables): k_cr = self.phase_param.k_cr(T) m_cr = self.phase_param.m_cr b_cr = self.phase_param.b_cr - stress_t_surf = variables[f"{Domain} {self.phase_name}particle surface tangential stress [Pa]"] + stress_t_surf = variables[ + f"{Domain} {self.phase_name}particle surface tangential stress [Pa]" + ] l_cr = variables[f"{Domain} {self.phase_name}particle crack length [m]"] # # compressive stress will not lead to crack propagation dK_SIF = stress_t_surf * b_cr * pybamm.sqrt(np.pi * l_cr) * (stress_t_surf >= 0) @@ -84,11 +86,17 @@ def set_rhs(self, variables): domain, Domain = self.domain_Domain if self.x_average is True: - l_cr = variables[f"X-averaged {domain} {self.phase_name}particle crack length [m]"] - dl_cr = variables[f"X-averaged {domain} {self.phase_name}particle cracking rate [m.s-1]"] + l_cr = variables[ + f"X-averaged {domain} {self.phase_name}particle crack length [m]" + ] + dl_cr = variables[ + f"X-averaged {domain} {self.phase_name}particle cracking rate [m.s-1]" + ] else: l_cr = variables[f"{Domain} {self.phase_name}particle crack length [m]"] - dl_cr = variables[f"{Domain} {self.phase_name}particle cracking rate [m.s-1]"] + dl_cr = variables[ + f"{Domain} {self.phase_name}particle cracking rate [m.s-1]" + ] self.rhs = {l_cr: dl_cr} def set_initial_conditions(self, variables): @@ -98,7 +106,9 @@ def set_initial_conditions(self, variables): if self.x_average is True: l_cr = variables[f"X-averaged {domain} particle crack length [m]"] else: - l_cr = variables[f"X-averaged {domain} {self.phase_name}particle crack length [m]"] + l_cr = variables[ + f"X-averaged {domain} {self.phase_name}particle crack length [m]" + ] l_cr_0 = pybamm.PrimaryBroadcast(l_cr_0, f"{domain} electrode") self.initial_conditions = {l_cr: l_cr_0} @@ -108,10 +118,12 @@ def add_events_from(self, variables): if self.x_average is True: l_cr = variables[f"X-averaged {domain} particle crack length [m]"] else: - l_cr = variables[f"X-averaged {domain} {self.phase_name}particle crack length [m]"] + l_cr = variables[ + f"X-averaged {domain} {self.phase_name}particle crack length [m]" + ] self.events.append( pybamm.Event( - f"{domain} {self.phase_name}particle crack length larger than particle radius", + f"{domain} {self.phase_name}particle crack length larger than particle radius", 1 - pybamm.max(l_cr) / self.phase_param.R_typ, ) ) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 1f297f6b1b..9b167cb535 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -37,7 +37,9 @@ def get_coupled_variables(self, variables): L_pl_k = variables[f"{Domain} lithium plating thickness [m]"] L_dead_k = variables[f"{Domain} dead lithium thickness [m]"] L_sei_cr_k = variables[f"{Domain} total SEI on cracks thickness [m]"] - roughness_k = variables[f"{Domain} electrode {self.phase_name}roughness ratio"] + roughness_k = variables[ + f"{Domain} electrode {self.phase_name}roughness ratio" + ] L_tot = ( (L_sei_k - L_sei_0) diff --git a/src/pybamm/parameters/lithium_ion_parameters.py b/src/pybamm/parameters/lithium_ion_parameters.py index 884e5ef1f3..e0c573c613 100644 --- a/src/pybamm/parameters/lithium_ion_parameters.py +++ b/src/pybamm/parameters/lithium_ion_parameters.py @@ -78,7 +78,6 @@ def _set_parameters(self): self.voltage_high_cut = self.elec.voltage_high_cut self.ocp_soc_0 = self.elec.ocp_soc_0 self.ocp_soc_100 = self.elec.ocp_soc_100 - # Domain parameters for domain in self.domain_params.values(): @@ -269,7 +268,6 @@ def _set_parameters(self): self.b_s = self.geo.b_s self.tau_s = self.geo.tau_s - # Utilisation parameters self.u_init = pybamm.Parameter( f"Initial {domain} electrode interface utilisation" @@ -294,7 +292,6 @@ def sigma(self, T): f"{Domain} electrode conductivity [S.m-1]", inputs ) - def LAM_rate_current(self, i, T): """ Dimensional rate of loss of active material as a function of applied current @@ -490,13 +487,17 @@ def _set_parameters(self): # Mechanical parameters self.nu = pybamm.Parameter(f"{pref}{Domain} electrode Poisson's ratio") self.c_0 = pybamm.Parameter( - f"{pref}{Domain} electrode reference concentration for free of deformation " - "[mol.m-3]" + f"{pref}{Domain} electrode reference concentration for free of deformation " + "[mol.m-3]" + ) + self.l_cr_0 = pybamm.Parameter( + f"{pref}{Domain} electrode initial crack length [m]" + ) + self.w_cr = pybamm.Parameter( + f"{pref}{Domain} electrode initial crack width [m]" ) - self.l_cr_0 = pybamm.Parameter(f"{pref}{Domain} electrode initial crack length [m]") - self.w_cr = pybamm.Parameter(f"{pref}{Domain} electrode initial crack width [m]") self.rho_cr = pybamm.Parameter( - f"{pref}{Domain} electrode number of cracks per unit area [m-2]" + f"{pref}{Domain} electrode number of cracks per unit area [m-2]" ) self.b_cr = pybamm.Parameter(f"{pref}{Domain} electrode Paris' law constant b") self.m_cr = pybamm.Parameter(f"{pref}{Domain} electrode Paris' law constant m") @@ -768,24 +769,34 @@ def dxdU(self, U, T): def Omega(self, sto, T): """Dimensional partial molar volume of Li in solid solution [m3.mol-1]""" Domain = self.domain.capitalize() - inputs = {f"{self.phase_prefactor}{Domain} particle stoichiometry": sto, "Temperature [K]": T} + inputs = { + f"{self.phase_prefactor}{Domain} particle stoichiometry": sto, + "Temperature [K]": T, + } return pybamm.FunctionParameter( - f"{self.phase_prefactor}{Domain} electrode partial molar volume [m3.mol-1]", inputs + f"{self.phase_prefactor}{Domain} electrode partial molar volume [m3.mol-1]", + inputs, ) + def E(self, sto, T): """Dimensional Young's modulus""" Domain = self.domain.capitalize() - inputs = {f"{self.phase_prefactor}{Domain} particle stoichiometry": sto, "Temperature [K]": T} + inputs = { + f"{self.phase_prefactor}{Domain} particle stoichiometry": sto, + "Temperature [K]": T, + } return pybamm.FunctionParameter( f"{self.phase_prefactor}{Domain} electrode Young's modulus [Pa]", inputs ) + def k_cr(self, T): """ Cracking rate for the electrode; """ Domain = self.domain.capitalize() return pybamm.FunctionParameter( - f"{self.phase_prefactor}{Domain} electrode cracking rate", {"Temperature [K]": T} + f"{self.phase_prefactor}{Domain} electrode cracking rate", + {"Temperature [K]": T}, ) def t_change(self, sto): @@ -799,5 +810,3 @@ def t_change(self, sto): f"{self.phase_prefactor}{Domain} electrode volume change", }, ) - - From dce8d03e3c65d47bb01d70ed49618f9371ebec6c Mon Sep 17 00:00:00 2001 From: mohammedasher Date: Wed, 23 Oct 2024 10:28:39 +0100 Subject: [PATCH 14/20] Fixed indentation and syntax errors in reaction_driven_porosity.py --- .../models/submodels/porosity/reaction_driven_porosity.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 8147f933d9..f172a1d782 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -42,14 +42,14 @@ def get_coupled_variables(self, variables): f"{Domain} electrode {self.phase_name}roughness ratio" ] - L_tot = ( + L_tot = ( (L_sei_k - L_sei_0) + L_pl_k + L_dead_k + L_sei_cr_k * (roughness_k - 1) ) - a_k = variables[ + a_k = variables[ f"{Domain} electrode {phase_name}" "surface area to volume ratio [m-1]" ] @@ -58,7 +58,7 @@ def get_coupled_variables(self, variables): # They could be included (e.g. for a sphere it is # a_n * (L_tot + L_tot ** 2 / R_n + L_tot ** # 3 / (3 * R_n ** 2))) # but it is not clear if it is relevant or not. - delta_eps_k += -a_k * L_tot + delta_eps_k += -a_k * L_tot domain_param = self.param.domain_params[domain.split()[0]] eps_k = domain_param.epsilon_init + delta_eps_k From 3996799ac85b71f9d2b7cde47e24d2f36886e740 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 23 Oct 2024 09:29:13 +0000 Subject: [PATCH 15/20] style: pre-commit fixes --- .../porosity/reaction_driven_porosity.py | 24 +++++++++---------- 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index f172a1d782..90418745dc 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -43,21 +43,21 @@ def get_coupled_variables(self, variables): ] L_tot = ( - (L_sei_k - L_sei_0) - + L_pl_k - + L_dead_k - + L_sei_cr_k * (roughness_k - 1) - ) + (L_sei_k - L_sei_0) + + L_pl_k + + L_dead_k + + L_sei_cr_k * (roughness_k - 1) + ) a_k = variables[ - f"{Domain} electrode {phase_name}" - "surface area to volume ratio [m-1]" - ] + f"{Domain} electrode {phase_name}" + "surface area to volume ratio [m-1]" + ] - # This assumes a thin film so curvature effects are neglected. - # They could be included (e.g. for a sphere it is - # a_n * (L_tot + L_tot ** 2 / R_n + L_tot ** # 3 / (3 * R_n ** 2))) - # but it is not clear if it is relevant or not. + # This assumes a thin film so curvature effects are neglected. + # They could be included (e.g. for a sphere it is + # a_n * (L_tot + L_tot ** 2 / R_n + L_tot ** # 3 / (3 * R_n ** 2))) + # but it is not clear if it is relevant or not. delta_eps_k += -a_k * L_tot domain_param = self.param.domain_params[domain.split()[0]] From 7576847ea1a2884ba86a4b87f2f40f67ce25f354 Mon Sep 17 00:00:00 2001 From: mohammedasher Date: Wed, 23 Oct 2024 10:32:30 +0100 Subject: [PATCH 16/20] phase_name to self.phase_name in reaction_driven_porosity.py --- .../models/submodels/porosity/reaction_driven_porosity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index f172a1d782..7d07f15650 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -50,7 +50,7 @@ def get_coupled_variables(self, variables): ) a_k = variables[ - f"{Domain} electrode {phase_name}" + f"{Domain} electrode {self.phase_name}" "surface area to volume ratio [m-1]" ] From 7b1832899d44d591941fb712d63070e84be7f250 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 23 Oct 2024 09:34:00 +0000 Subject: [PATCH 17/20] style: pre-commit fixes --- .../models/submodels/porosity/reaction_driven_porosity.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 846bc8a4ce..0df9c830a7 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -50,9 +50,9 @@ def get_coupled_variables(self, variables): ) a_k = variables[ - f"{Domain} electrode {self.phase_name}" - "surface area to volume ratio [m-1]" - ] + f"{Domain} electrode {self.phase_name}" + "surface area to volume ratio [m-1]" + ] # This assumes a thin film so curvature effects are neglected. # They could be included (e.g. for a sphere it is From 07ac6fabb664d0e30f5a01e4c1b4dd76de7dd324 Mon Sep 17 00:00:00 2001 From: mohammedasher Date: Wed, 23 Oct 2024 10:40:33 +0100 Subject: [PATCH 18/20] Applied pre-commit fixes after removing ghost files --- .../models/submodels/porosity/reaction_driven_porosity.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py index 846bc8a4ce..0df9c830a7 100644 --- a/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py +++ b/src/pybamm/models/submodels/porosity/reaction_driven_porosity.py @@ -50,9 +50,9 @@ def get_coupled_variables(self, variables): ) a_k = variables[ - f"{Domain} electrode {self.phase_name}" - "surface area to volume ratio [m-1]" - ] + f"{Domain} electrode {self.phase_name}" + "surface area to volume ratio [m-1]" + ] # This assumes a thin film so curvature effects are neglected. # They could be included (e.g. for a sphere it is From d8d2b6c78238966555201940f509d20b71de9cf9 Mon Sep 17 00:00:00 2001 From: "Mohammed Asheruddin (Asher)" <168521559+mohammedasher@users.noreply.github.com> Date: Fri, 15 Nov 2024 15:22:22 +0000 Subject: [PATCH 19/20] Update lithium_ion_parameters.py --- src/pybamm/parameters/lithium_ion_parameters.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/pybamm/parameters/lithium_ion_parameters.py b/src/pybamm/parameters/lithium_ion_parameters.py index e0c573c613..b332edc1bc 100644 --- a/src/pybamm/parameters/lithium_ion_parameters.py +++ b/src/pybamm/parameters/lithium_ion_parameters.py @@ -806,7 +806,5 @@ def t_change(self, sto): Domain = self.domain.capitalize() return pybamm.FunctionParameter( f"{self.phase_prefactor}{Domain} electrode volume change", - { - f"{self.phase_prefactor}{Domain} electrode volume change", - }, + {f"{self.phase_prefactor}{Domain} electrode stoichiometry": sto}, ) From 3adba7e1764f46042744e7362757dcfd889011fb Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 18 Dec 2024 17:27:28 +0000 Subject: [PATCH 20/20] style: pre-commit fixes --- src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py | 2 -- src/pybamm/models/submodels/interface/sei/base_sei.py | 2 -- 2 files changed, 4 deletions(-) diff --git a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py index 124daa1792..ff5722f7e2 100644 --- a/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py +++ b/src/pybamm/input/parameters/lithium_ion/Chen2020_composite.py @@ -535,7 +535,6 @@ def get_parameter_values(): return { "chemistry": "lithium_ion", # sei - "Primary: Negative ratio of lithium moles to SEI moles": 2.0, "Primary: Negative inner SEI reaction proportion": 0.5, "Primary: Negative inner SEI partial molar volume [m3.mol-1]": 9.585e-05, @@ -577,7 +576,6 @@ def get_parameter_values(): "Secondary: Negative SEI open-circuit potential [V]": 0.4, "Secondary: Negative SEI growth activation energy [J.mol-1]": 0.0, "Negative electrode reaction-driven LAM factor [m3.mol-1]": 0.0, - "Positive electrode reaction-driven LAM factor [m3.mol-1]": 0.0, # cell "Negative current collector thickness [m]": 1.2e-05, diff --git a/src/pybamm/models/submodels/interface/sei/base_sei.py b/src/pybamm/models/submodels/interface/sei/base_sei.py index 8f1d95f42a..2d88fba60e 100644 --- a/src/pybamm/models/submodels/interface/sei/base_sei.py +++ b/src/pybamm/models/submodels/interface/sei/base_sei.py @@ -184,14 +184,12 @@ def _get_standard_concentration_variables(self, variables): ) # Concentration variables are handled slightly differently for SEI on cracks elif self.reaction == "SEI on cracks": - L_inner_cr = variables[f"{Domain} inner {reaction_name}thickness [m]"] L_outer_cr = variables[f"{Domain} outer {reaction_name}thickness [m]"] roughness = variables[ f"{Domain} electrode {self.phase_name}roughness ratio" ] - n_SEI_cr = L_cr * L_to_n * (roughness - 1) # SEI on cracks concentration n_SEI_cr_xav = pybamm.x_average(n_SEI_cr) n_SEI_cr_av = pybamm.yz_average(n_SEI_cr_xav)