update pouch models

pybamm/models/submodels/current_collector/potential_pair.py

@@ -80,12 +80,13 @@ def set_boundary_conditions(self, variables):
 
         param = self.param
         applied_current_density = variables["Total current density [A.m-2]"]
-        cc_area = self._get_effective_current_collector_area()
+        total_current = applied_current_density * param.A_cc
 
-        # cc_area appears here due to choice of non-dimensionalisation
-        pos_tab_bc = (
-            -applied_current_density * cc_area / (param.p.sigma_cc * param.p.L_cc)
-        )
+        # In the 1+1D model, the behaviour is averaged over the y-direction, so the
+        # effective tab area is the cell width multiplied by the current collector
+        # thickness
+        positive_tab_area = param.L_y * param.p.L_cc
+        pos_tab_bc = -total_current / (param.p.sigma_cc * positive_tab_area)
 
         # Boundary condition needs to be on the variables that go into the Laplacian,
         # even though phi_s_cp isn't a pybamm.Variable object
@@ -100,10 +101,6 @@ def set_boundary_conditions(self, variables):
             },
         }
 
-    def _get_effective_current_collector_area(self):
-        """In the 1+1D models the current collector effectively has surface area l_z"""
-        return self.param.L_z
-
 
 class PotentialPair2plus1D(BasePotentialPair):
     """Base class for a 2+1D potential pair model"""
@@ -117,21 +114,18 @@ def set_boundary_conditions(self, variables):
 
         param = self.param
         applied_current_density = variables["Total current density [A.m-2]"]
-        cc_area = self._get_effective_current_collector_area()
+        total_current = applied_current_density * param.A_cc
 
         # Note: we divide by the *numerical* tab area so that the correct total
         # current is applied. That is, numerically integrating the current density
         # around the boundary gives the applied current exactly.
-
         positive_tab_area = pybamm.BoundaryIntegral(
             pybamm.PrimaryBroadcast(param.p.L_cc, "current collector"),
             region="positive tab",
         )
 
         # cc_area appears here due to choice of non-dimensionalisation
-        pos_tab_bc = (
-            -applied_current_density * cc_area / (param.p.sigma_cc * positive_tab_area)
-        )
+        pos_tab_bc = -total_current / (param.p.sigma_cc * positive_tab_area)
 
         # Boundary condition needs to be on the variables that go into the Laplacian,
         # even though phi_s_cp isn't a pybamm.Variable object
@@ -160,7 +154,3 @@ def set_boundary_conditions(self, variables):
                 "positive tab": (pos_tab_bc, "Neumann"),
             },
         }
-
-    def _get_effective_current_collector_area(self):
-        """Return the area of the current collector."""
-        return self.param.L_y * self.param.L_z

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

@@ -82,7 +82,7 @@ def set_rhs(self, variables):
 
         self.rhs = {
             T_av: (
-                pybamm.laplacian(T_av)
+                pybamm.div(self.param.lambda_eff(T_av) * pybamm.grad(T_av))
                 + Q_av
                 + total_cooling_coefficient * (T_av - T_amb)
             )

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

@@ -82,9 +82,11 @@ def set_rhs(self, variables):
         # correct mass matrix when discretised. The first argument is the source term
         # and the second argument is the variable governed by the equation that the
         # source term appears in.
+        # Note: not correct if lambda_eff is a function of T_av - need to implement div
+        # in 2D rather than doing laplacian directly
         self.rhs = {
             T_av: (
-                pybamm.laplacian(T_av)
+                self.param.lambda_eff(T_av) * pybamm.laplacian(T_av)
                 + pybamm.source(Q_av, T_av)
                 + pybamm.source(yz_surface_cooling_coefficient * (T_av - T_amb), T_av)
                 + pybamm.source(

pybamm/parameters/lithium_ion_parameters.py

@@ -63,6 +63,7 @@ def _set_parameters(self):
         self.h_edge = self.therm.h_edge
         self.h_total = self.therm.h_total
         self.rho_c_p_eff = self.therm.rho_c_p_eff
+        self.lambda_eff = self.therm.lambda_eff
 
         # Macroscale geometry
         self.L_x = self.geo.L_x

pybamm/parameters/thermal_parameters.py

@@ -72,6 +72,16 @@ def rho_c_p_eff(self, T):
             + self.p.rho_c_p_cc(T) * self.geo.p.L_cc
         ) / self.geo.L
 
+    def lambda_eff(self, T):
+        """Effective thermal conductivity [W.m-1.K-1]"""
+        return (
+            self.n.lambda_cc(T) * self.geo.n.L_cc
+            + self.n.lambda_(T) * self.geo.n.L
+            + self.s.lambda_(T) * self.geo.s.L
+            + self.p.lambda_(T) * self.geo.p.L
+            + self.p.lambda_cc(T) * self.geo.p.L_cc
+        ) / self.geo.L
+
 
 class DomainThermalParameters(BaseParameters):
     def __init__(self, domain, main_param):