Issue 1571 1 d cc opposite tabs

CHANGELOG.md

@@ -2,6 +2,7 @@
 
 ## Features
 
+-   Tabs can now be placed at the bottom of the cell in 1+1D thermal models ([#1581](https://github.com/pybamm-team/PyBaMM/pull/1581)) 
 -   Added temperature dependence on electrode electronic conductivity ([#1570](https://github.com/pybamm-team/PyBaMM/pull/1570))
 -   Added a new lithium-ion model `MPM` or Many-Particle Model, with a distribution of particle sizes in each electrode. ([#1529](https://github.com/pybamm-team/PyBaMM/pull/1529))
 -   Added 2 new submodels for lithium transport in a size distribution of electrode particles: Fickian diffusion (`FickianSingleSizeDistribution`) and uniform concentration profile (`FastSingleSizeDistribution`). ([#1529](https://github.com/pybamm-team/PyBaMM/pull/1529))

pybamm/expression_tree/binary_operators.py

@@ -502,6 +502,38 @@ def inner(left, right):
     return pybamm.simplify_if_constant(pybamm.Inner(left, right))
 
 
+class Equality(BinaryOperator):
+    """
+    A node in the expression tree representing an equality comparison between two
+    nodes. Returns 1 if the two nodes evaluate to the same thing and 0 otherwise.
+    **Extends:** :class:`BinaryOperator`
+    """
+
+    def __init__(self, left, right):
+        """See :meth:`pybamm.BinaryOperator.__init__()`."""
+        super().__init__("==", left, right)
+
+    def diff(self, variable):
+        """See :meth:`pybamm.Symbol.diff()`."""
+        # Equality should always be multiplied by something else so hopefully don't
+        # need to worry about shape
+        return pybamm.Scalar(0)
+
+    def _binary_jac(self, left_jac, right_jac):
+        """See :meth:`pybamm.BinaryOperator._binary_jac()`."""
+        # Equality should always be multiplied by something else so hopefully don't
+        # need to worry about shape
+        return pybamm.Scalar(0)
+
+    def _binary_evaluate(self, left, right):
+        """See :meth:`pybamm.BinaryOperator._binary_evaluate()`."""
+        return int(left == right)
+
+    def _binary_new_copy(self, left, right):
+        """See :meth:`pybamm.BinaryOperator._binary_new_copy()`."""
+        return pybamm.Equality(left, right)
+
+
 class _Heaviside(BinaryOperator):
     """
     A node in the expression tree representing a heaviside step function.

pybamm/models/full_battery_models/base_battery_model.py

@@ -3,7 +3,6 @@
 #
 
 import pybamm
-import warnings
 
 
 class BatteryModelOptions(pybamm.FuzzyDict):
@@ -332,12 +331,6 @@ def __init__(self, extra_options):
                 "Use 'uniform profile' instead."
             )
 
-        if options["thermal"] == "x-lumped" and options["dimensionality"] == 1:
-            warnings.warn(
-                "1+1D Thermal models are only valid if both tabs are "
-                "placed at the top of the cell."
-            )
-
         for option, value in options.items():
             if option == "external submodels" or option == "working electrode":
                 pass

pybamm/models/submodels/particle/size_distribution/fast_single_distribution.py

@@ -95,12 +95,11 @@ def get_fundamental_variables(self):
         # quantities. Necessary for output variables "Total lithium in
         # negative electrode [mol]", etc, to be calculated correctly
         f_v_dist = variables[
-            "X-averaged " + self.domain.lower()
+            "X-averaged "
+            + self.domain.lower()
             + " volume-weighted particle-size distribution"
         ]
-        c_s_surf_xav = pybamm.Integral(
-            f_v_dist * c_s_surf_xav_distribution, R
-        )
+        c_s_surf_xav = pybamm.Integral(f_v_dist * c_s_surf_xav_distribution, R)
         c_s_xav = pybamm.PrimaryBroadcast(
             c_s_surf_xav, [self.domain.lower() + " particle"]
         )
@@ -164,27 +163,3 @@ def set_initial_conditions(self, variables):
             c_init = self.param.c_p_init(1)
 
         self.initial_conditions = {c_s_surf_xav_distribution: c_init}
-
-    def set_events(self, variables):
-        c_s_surf_xav_distribution = variables[
-            "X-averaged "
-            + self.domain.lower()
-            + " particle surface concentration distribution"
-        ]
-        tol = 1e-4
-
-        self.events.append(
-            pybamm.Event(
-                "Minimum " + self.domain.lower() + " particle surface concentration",
-                pybamm.min(c_s_surf_xav_distribution) - tol,
-                pybamm.EventType.TERMINATION,
-            )
-        )
-
-        self.events.append(
-            pybamm.Event(
-                "Maximum " + self.domain.lower() + " particle surface concentration",
-                (1 - tol) - pybamm.max(c_s_surf_xav_distribution),
-                pybamm.EventType.TERMINATION,
-            )
-        )

pybamm/models/submodels/particle/size_distribution/fickian_single_distribution.py

@@ -84,14 +84,13 @@ def get_fundamental_variables(self):
         # quantities. Necessary for output variables "Total lithium in
         # negative electrode [mol]", etc, to be calculated correctly
         f_v_dist = variables[
-            "X-averaged " + self.domain.lower()
+            "X-averaged "
+            + self.domain.lower()
             + " volume-weighted particle-size distribution"
         ]
         c_s_xav = pybamm.Integral(f_v_dist * c_s_xav_distribution, R)
         c_s = pybamm.SecondaryBroadcast(c_s_xav, [self.domain.lower() + " electrode"])
-        variables.update(
-            self._get_standard_concentration_variables(c_s, c_s_xav)
-        )
+        variables.update(self._get_standard_concentration_variables(c_s, c_s_xav))
         return variables
 
     def get_coupled_variables(self, variables):
@@ -105,7 +104,10 @@ def get_coupled_variables(self, variables):
             variables["X-averaged " + self.domain.lower() + " electrode temperature"],
             [self.domain.lower() + " particle size"],
         )
-        T_k_xav = pybamm.PrimaryBroadcast(T_k_xav, [self.domain.lower() + " particle"],)
+        T_k_xav = pybamm.PrimaryBroadcast(
+            T_k_xav,
+            [self.domain.lower() + " particle"],
+        )
 
         if self.domain == "Negative":
             N_s_xav_distribution = -self.param.D_n(
@@ -119,11 +121,13 @@ def get_coupled_variables(self, variables):
         # Standard R-averaged flux variables. Average using the area-weighted
         # distribution
         f_a_dist = variables[
-            "X-averaged " + self.domain.lower()
+            "X-averaged "
+            + self.domain.lower()
             + " area-weighted particle-size distribution"
         ]
         f_a_dist = pybamm.PrimaryBroadcast(
-            f_a_dist, [self.domain.lower() + " particle"],
+            f_a_dist,
+            [self.domain.lower() + " particle"],
         )
         # must use "R_spatial_variable" as integration variable, since "R" is a
         # broadcast
@@ -157,7 +161,8 @@ def set_rhs(self, variables):
         # Spatial variable R, broadcast into particle
         R_spatial_variable = variables[self.domain + " particle sizes"]
         R = pybamm.PrimaryBroadcast(
-            R_spatial_variable, [self.domain.lower() + " particle"],
+            R_spatial_variable,
+            [self.domain.lower() + " particle"],
         )
         if self.domain == "Negative":
             self.rhs = {
@@ -243,28 +248,3 @@ def set_initial_conditions(self, variables):
             c_init = self.param.c_p_init(1)
 
         self.initial_conditions = {c_s_xav_distribution: c_init}
-
-    def set_events(self, variables):
-        c_s_surf_xav_distribution = variables[
-            "X-averaged "
-            + self.domain.lower()
-            + " particle surface concentration distribution"
-        ]
-        tol = 1e-4
-
-        self.events.append(
-            pybamm.Event(
-                "Minimum " + self.domain.lower() + " particle surface concentration",
-                pybamm.min(c_s_surf_xav_distribution) - tol,
-                pybamm.EventType.TERMINATION,
-            )
-        )
-
-        self.events.append(
-            pybamm.Event(
-                "Maximum " + self.domain.lower() + " particle surface concentration",
-                (1 - tol) - pybamm.max(c_s_surf_xav_distribution),
-                pybamm.EventType.TERMINATION,
-            )
-        )
-

pybamm/models/submodels/thermal/pouch_cell/pouch_cell_1D_current_collectors.py

@@ -90,51 +90,68 @@ def set_rhs(self, variables):
         }
 
     def set_boundary_conditions(self, variables):
+        param = self.param
         T_amb = variables["Ambient temperature"]
         T_av = variables["X-averaged cell temperature"]
         T_av_top = pybamm.boundary_value(T_av, "right")
         T_av_bottom = pybamm.boundary_value(T_av, "left")
 
-        # Tab cooling only implemented for both tabs at the top.
-        negative_tab_area = self.param.l_tab_n * self.param.l_cn
-        positive_tab_area = self.param.l_tab_p * self.param.l_cp
-        total_top_area = self.param.l * self.param.l_y
-        non_tab_top_area = total_top_area - negative_tab_area - positive_tab_area
-
-        negative_tab_cooling_coefficient = (
-            self.param.h_tab_n / self.param.delta * negative_tab_area / total_top_area
-        )
-        positive_tab_cooling_coefficient = (
-            self.param.h_tab_p / self.param.delta * positive_tab_area / total_top_area
+        # find tab locations (top vs bottom)
+        l_z = param.l_z
+        neg_tab_z = param.centre_z_tab_n
+        pos_tab_z = param.centre_z_tab_p
+        neg_tab_top_bool = pybamm.Equality(neg_tab_z, l_z)
+        neg_tab_bottom_bool = pybamm.Equality(neg_tab_z, 0)
+        pos_tab_top_bool = pybamm.Equality(pos_tab_z, l_z)
+        pos_tab_bottom_bool = pybamm.Equality(pos_tab_z, 0)
+
+        # calculate tab vs non-tab area on top and bottom
+        neg_tab_area = param.l_tab_n * param.l_cn
+        pos_tab_area = param.l_tab_p * param.l_cp
+        total_area = param.l * param.l_y
+
+        non_tab_top_area = (
+            total_area
+            - neg_tab_area * neg_tab_top_bool
+            - pos_tab_area * pos_tab_top_bool
         )
-
-        top_edge_cooling_coefficient = (
-            self.param.h_edge / self.param.delta * non_tab_top_area / total_top_area
+        non_tab_bottom_area = (
+            total_area
+            - neg_tab_area * neg_tab_bottom_bool
+            - pos_tab_area * pos_tab_bottom_bool
         )
 
-        bottom_edge_cooling_coefficient = (
-            self.param.h_edge / self.param.delta * total_top_area / total_top_area
+        # calculate effective cooling coefficients
+        top_cooling_coefficient = (
+            (
+                param.h_tab_n * neg_tab_area * neg_tab_top_bool
+                + param.h_tab_p * pos_tab_area * pos_tab_top_bool
+                + param.h_edge * non_tab_top_area
+            )
+            / param.delta
+            / total_area
         )
-
-        total_top_cooling_coefficient = (
-            negative_tab_cooling_coefficient
-            + positive_tab_cooling_coefficient
-            + top_edge_cooling_coefficient
+        bottom_cooling_coefficient = (
+            (
+                param.h_tab_n * neg_tab_area * neg_tab_bottom_bool
+                + param.h_tab_p * pos_tab_area * pos_tab_bottom_bool
+                + param.h_edge * non_tab_bottom_area
+            )
+            / param.delta
+            / total_area
         )
 
-        total_bottom_cooling_coefficient = bottom_edge_cooling_coefficient
-
         # just use left and right for clarity
         # left = bottom of cell (z=0)
         # right = top of cell (z=L_z)
         self.boundary_conditions = {
             T_av: {
                 "left": (
-                    total_bottom_cooling_coefficient * (T_av_bottom - T_amb),
+                    bottom_cooling_coefficient * (T_av_bottom - T_amb),
                     "Neumann",
                 ),
                 "right": (
-                    -total_top_cooling_coefficient * (T_av_top - T_amb),
+                    -top_cooling_coefficient * (T_av_top - T_amb),
                     "Neumann",
                 ),
             }

tests/unit/test_expression_tree/test_binary_operators.py

@@ -323,6 +323,14 @@ def test_heaviside(self):
         self.assertEqual(heav.evaluate(y=np.array([0])), 1)
         self.assertEqual(str(heav), "y[0:1] <= 1.0")
 
+    def test_equality(self):
+        a = pybamm.Scalar(1)
+        b = pybamm.StateVector(slice(0, 1))
+        equal = pybamm.Equality(a, b)
+        self.assertEqual(equal.evaluate(y=np.array([1])), 1)
+        self.assertEqual(equal.evaluate(y=np.array([2])), 0)
+        self.assertEqual(str(equal), "1.0 == y[0:1]")
+
     def test_sigmoid(self):
         a = pybamm.Scalar(1)
         b = pybamm.StateVector(slice(0, 1))

tests/unit/test_expression_tree/test_operations/test_copy.py

@@ -58,6 +58,7 @@ def test_symbol_new_copy(self):
             pybamm.minimum(a, b),
             pybamm.maximum(a, b),
             pybamm.SparseStack(mat, mat),
+            pybamm.Equality(a, b),
         ]:
             self.assertEqual(symbol.id, symbol.new_copy().id)
 

tests/unit/test_expression_tree/test_operations/test_jac.py

@@ -306,6 +306,14 @@ def test_jac_of_heaviside(self):
             ((a < y) * y ** 2).jac(y).evaluate(y=-5 * np.ones(5)).toarray(), 0
         )
 
+    def test_jac_of_equality(self):
+        a = pybamm.Scalar(1)
+        y = pybamm.StateVector(slice(0, 1))
+        np.testing.assert_array_equal(
+            (pybamm.Equality(a, y)).jac(y).evaluate(),
+            0,
+        )
+
     def test_jac_of_modulo(self):
         a = pybamm.Scalar(3)
         y = pybamm.StateVector(slice(0, 5))