Add half cell submodels

docs/source/models/submodels/electrode/ohm/index.rst

@@ -8,3 +8,4 @@ Ohmic
   composite_ohm
   full_ohm
   surface_form_ohm
+  li_metal_explicit

docs/source/models/submodels/electrode/ohm/li_metal_explicit.rst

@@ -0,0 +1,5 @@
+Explicit potential drop for lithium metal
+=========================================
+
+.. autoclass:: pybamm.electrode.ohm.LithiumMetalExplicit
+    :members:

examples/scripts/DFN.py

@@ -8,7 +8,7 @@
 pybamm.set_logging_level("INFO")
 
 # load model
-model = pybamm.lithium_ion.DFN()
+model = pybamm.lithium_ion.DFN({"working electrode": "positive"})
 # create geometry
 geometry = model.default_geometry
 

examples/scripts/DFN_half_cell.py

@@ -9,58 +9,51 @@
 
 # load model
 options = {"working electrode": "positive"}
-model = pybamm.lithium_ion.BasicDFNHalfCell(options=options)
+model1 = pybamm.lithium_ion.DFN(options=options)
+model2 = pybamm.lithium_ion.BasicDFNHalfCell(options=options)
 
-# create geometry
-geometry = model.default_geometry
+sols = []
+for model in [model1, model2]:
+    # for model in [model2]:
+    # create geometry
+    geometry = model.default_geometry
 
-# load parameter values
-chemistry = pybamm.parameter_sets.Chen2020
-param = pybamm.ParameterValues(chemistry=chemistry)
+    # load parameter values
+    chemistry = pybamm.parameter_sets.Xu2019
+    param = pybamm.ParameterValues(chemistry=chemistry)
 
-# add lithium counter electrode parameter values
-param.update(
-    {
-        "Lithium counter electrode exchange-current density [A.m-2]": 12.6,
-        "Lithium counter electrode conductivity [S.m-1]": 1.0776e7,
-        "Lithium counter electrode thickness [m]": 250e-6,
-    },
-    check_already_exists=False,
-)
-
-param["Initial concentration in negative electrode [mol.m-3]"] = 1000
-param["Current function [A]"] = 2.5
-
-# process model and geometry
-param.process_model(model)
-param.process_geometry(geometry)
+    # process model and geometry
+    param.process_model(model)
+    param.process_geometry(geometry)
 
-# set mesh
-var_pts = model.default_var_pts
-mesh = pybamm.Mesh(geometry, model.default_submesh_types, var_pts)
+    # set mesh
+    var_pts = model.default_var_pts
+    mesh = pybamm.Mesh(geometry, model.default_submesh_types, var_pts)
 
-# discretise model
-disc = pybamm.Discretisation(mesh, model.default_spatial_methods)
-disc.process_model(model)
+    # discretise model
+    disc = pybamm.Discretisation(mesh, model.default_spatial_methods)
+    disc.process_model(model)
 
-# solve model
-t_eval = np.linspace(0, 7200, 1000)
-solver = pybamm.CasadiSolver(mode="safe", atol=1e-6, rtol=1e-3)
-solution = solver.solve(model, t_eval)
+    # solve model
+    t_eval = np.linspace(0, 7200, 1000)
+    solver = pybamm.CasadiSolver(mode="safe", atol=1e-6, rtol=1e-3)
+    solution = solver.solve(model, t_eval)
+    sols.append(solution)
 
 # plot
 plot = pybamm.QuickPlot(
-    solution,
+    sols,
     [
-        "Working particle concentration [mol.m-3]",
+        # "Positive particle concentration [mol.m-3]",
         "Electrolyte concentration [mol.m-3]",
         "Current [A]",
-        "Working electrode potential [V]",
+        "Positive electrode potential [V]",
         "Electrolyte potential [V]",
-        "Total electrolyte concentration",
-        "Total lithium in working electrode [mol]",
-        "Working electrode open circuit potential [V]",
-        ["Terminal voltage [V]", "Voltage drop in the cell [V]"],
+        "Total lithium in electrolyte [mol]",
+        # "Total lithium in positive electrode [mol]",
+        "Positive electrode open circuit potential [V]",
+        ["Terminal voltage [V]"],  # , "Voltage drop in the cell [V]"],
+        "Negative electrode potential drop [V]",
     ],
     time_unit="seconds",
     spatial_unit="um",

pybamm/expression_tree/concatenations.py

@@ -397,6 +397,8 @@ def intersect(s1, s2):
 
 def simplified_concatenation(*children):
     """Perform simplifications on a concatenation."""
+    # remove children that are None
+    children = list(filter(lambda x: x is not None, children))
     # Create Concatenation to easily read domains
     concat = Concatenation(*children)
     # Simplify concatenation of broadcasts all with the same child to a single

pybamm/expression_tree/unary_operators.py

@@ -1265,17 +1265,19 @@ def x_average(symbol):
     if isinstance(symbol, (pybamm.PrimaryBroadcast, pybamm.FullBroadcast)):
         return symbol.reduce_one_dimension()
     # If symbol is a concatenation of Broadcasts, its average value is its child
-    elif (
-        isinstance(symbol, pybamm.Concatenation)
-        and all(isinstance(child, pybamm.Broadcast) for child in symbol.children)
-        and symbol.domain == ["negative electrode", "separator", "positive electrode"]
+    elif isinstance(symbol, pybamm.Concatenation) and all(
+        isinstance(child, pybamm.Broadcast) for child in symbol.children
     ):
-        a, b, c = [orp.orphans[0] for orp in symbol.orphans]
         geo = pybamm.geometric_parameters
         l_n = geo.l_n
         l_s = geo.l_s
         l_p = geo.l_p
-        out = (l_n * a + l_s * b + l_p * c) / (l_n + l_s + l_p)
+        if symbol.domain == ["negative electrode", "separator", "positive electrode"]:
+            a, b, c = [orp.orphans[0] for orp in symbol.orphans]
+            out = (l_n * a + l_s * b + l_p * c) / (l_n + l_s + l_p)
+        elif symbol.domain == ["separator", "positive electrode"]:
+            b, c = [orp.orphans[0] for orp in symbol.orphans]
+            out = (l_s * b + l_p * c) / (l_s + l_p)
         # To respect domains we may need to broadcast the child back out
         child = symbol.children[0]
         # If symbol being returned doesn't have empty domain, return it

pybamm/input/parameters/lithium_ion/cells/li_metal_Xu2019/parameters.csv

@@ -9,6 +9,6 @@ Electrode height [m],0.01,,guess to get correct area
 Electrode width [m],0.0154,,guess to get correct area
 ,,,
 # Electrical,,,
-Nominal cell capacity [A.h],0.024,,C/4 is I=0.6 mA
-Typical current [A],0.024,,1C current
-Current function [A],0.024,,1C current
+Nominal cell capacity [A.h],0.0024,,C/4 is I=0.6 mA
+Typical current [A],0.0024,,1C current
+Current function [A],0.0024,,1C current

pybamm/input/parameters/lithium_ion/experiments/1C_discharge_from_full_Xu2019/parameters.csv

@@ -2,6 +2,7 @@ Name [units],Value,Reference,Notes
 # Empty rows and rows starting with ‘#’ will be ignored,,,
 ,,,
 # Temperature
+Ambient temperature [K],298.15,25C,
 Reference temperature [K],298.15,25C,
 Heat transfer coefficient [W.m-2.K-1],10,,
 ,,,

pybamm/input/parameters/lithium_ion/negative_electrodes/li_metal_Xu2019/parameters.csv

@@ -3,14 +3,15 @@ Name [units],Value,Reference,Notes
 ,,,
 # Electrode properties,,,
 Negative electrode OCP [V],0,lithium metal,
+Negative electrode conductivity [S.m-1],1.0776e7,li metal conductivity not provided by Xu2019
 ,,,
 # Thermal parameters,,,
 Negative electrode OCP entropic change [V.K-1],0,,
 ,,,
 # Interfacial reactions,,,
 Negative electrode cation signed stoichiometry,-1,,
 Negative electrode electrons in reaction,1,,
-Lithium metal concentration [mol.m-3],76900,,
-Negative electrode exchange-current density [A.m-2],[function]li_metal_electrolyte_exchange_current_density_Xu2019,,
+Typical plated lithium concentration [mol.m-3],76900,,
+Exchange-current density for plating [A.m-2],[function]li_metal_electrolyte_exchange_current_density_Xu2019,,
 Negative electrode charge transfer coefficient,0.5,,
 Negative electrode double-layer capacity [F.m-2],0.2,,

pybamm/models/full_battery_models/base_battery_model.py

@@ -765,16 +765,18 @@ def set_external_circuit_submodel(self):
 
     def set_tortuosity_submodels(self):
         self.submodels["electrolyte tortuosity"] = pybamm.tortuosity.Bruggeman(
-            self.param, "Electrolyte"
+            self.param, "Electrolyte", self.options
         )
         self.submodels["electrode tortuosity"] = pybamm.tortuosity.Bruggeman(
-            self.param, "Electrode"
+            self.param, "Electrode", self.options
         )
 
     def set_thermal_submodel(self):
 
         if self.options["thermal"] == "isothermal":
-            thermal_submodel = pybamm.thermal.isothermal.Isothermal(self.param)
+            thermal_submodel = pybamm.thermal.isothermal.Isothermal(
+                self.param, self.options
+            )
 
         elif self.options["thermal"] == "lumped":
             thermal_submodel = pybamm.thermal.Lumped(
@@ -940,9 +942,7 @@ def set_voltage_variables(self):
         )
         # Variables for calculating the equivalent circuit model (ECM) resistance
         # Need to compare OCV to initial value to capture this as an overpotential
-        ocv_init = self.param.U_p(
-            self.param.c_p_init(1), self.param.T_init
-        ) - self.param.U_n(self.param.c_n_init(0), self.param.T_init)
+        ocv_init = self.param.U_p_init - self.param.U_n_init
         ocv_init_dim = (
             self.param.U_p_ref
             - self.param.U_n_ref