Merge pull request #2758 from pybamm-team/issue-2666-electrolyte_conservation

#2666 put t_+ back inside div, loosen tolerance for test

Author: valentinsulzer

CHANGELOG.md

@@ -9,6 +9,7 @@
 
 ## Bug fixes
 
+- Fixed electrolyte conservation in the case of concentration-dependent transference number ([#2758](https://github.com/pybamm-team/PyBaMM/pull/2758))
 - Fixed `plot_voltage_components` so that the sum of overpotentials is now equal to the voltage ([#2740](https://github.com/pybamm-team/PyBaMM/pull/2740))
 
 ## Optimizations

examples/notebooks/Getting Started/Tutorial 7 - Model options.ipynb

@@ -59,7 +59,7 @@
     {
      "data": {
       "text/plain": [
-       "<pybamm.solvers.solution.Solution at 0x7f29f1333e50>"
+       "<pybamm.solvers.solution.Solution at 0x139461520>"
       ]
      },
      "execution_count": 3,
@@ -70,10 +70,7 @@
    "source": [
     "model = pybamm.lithium_ion.SPMe(options=options) # loading in options\n",
     "\n",
-    "parameter_values = model.default_parameter_values\n",
-    "parameter_values[\"Current function [A]\"] = 3\n",
-    "\n",
-    "sim = pybamm.Simulation(model, parameter_values=parameter_values)\n",
+    "sim = pybamm.Simulation(model)\n",
     "sim.solve([0, 3600])"
    ]
   },
@@ -92,16 +89,26 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "c00d126a3fc543829596dd24f8b9d626",
+       "model_id": "0fc5db759e804d3fb2793222ddbec5d4",
        "version_major": 2,
        "version_minor": 0
       },
       "text/plain": [
-       "interactive(children=(FloatSlider(value=0.0, description='t', max=800.0000000000001, step=8.000000000000002), …"
+       "interactive(children=(FloatSlider(value=0.0, description='t', max=1.0, step=0.01), Output()), _dom_classes=('w…"
       ]
      },
      "metadata": {},
      "output_type": "display_data"
+    },
+    {
+     "data": {
+      "text/plain": [
+       "<pybamm.plotting.quick_plot.QuickPlot at 0x139461fd0>"
+      ]
+     },
+     "execution_count": 4,
+     "metadata": {},
+     "output_type": "execute_result"
     }
    ],
    "source": [
@@ -139,7 +146,8 @@
       "[1] Joel A. E. Andersson, Joris Gillis, Greg Horn, James B. Rawlings, and Moritz Diehl. CasADi – A software framework for nonlinear optimization and optimal control. Mathematical Programming Computation, 11(1):1–36, 2019. doi:10.1007/s12532-018-0139-4.\n",
       "[2] Charles R. Harris, K. Jarrod Millman, Stéfan J. van der Walt, Ralf Gommers, Pauli Virtanen, David Cournapeau, Eric Wieser, Julian Taylor, Sebastian Berg, Nathaniel J. Smith, and others. Array programming with NumPy. Nature, 585(7825):357–362, 2020. doi:10.1038/s41586-020-2649-2.\n",
       "[3] Scott G. Marquis, Valentin Sulzer, Robert Timms, Colin P. Please, and S. Jon Chapman. An asymptotic derivation of a single particle model with electrolyte. Journal of The Electrochemical Society, 166(15):A3693–A3706, 2019. doi:10.1149/2.0341915jes.\n",
-      "[4] Valentin Sulzer, Scott G. Marquis, Robert Timms, Martin Robinson, and S. Jon Chapman. Python Battery Mathematical Modelling (PyBaMM). ECSarXiv. February, 2020. doi:10.1149/osf.io/67ckj.\n",
+      "[4] Valentin Sulzer, Scott G. Marquis, Robert Timms, Martin Robinson, and S. Jon Chapman. Python Battery Mathematical Modelling (PyBaMM). Journal of Open Research Software, 9(1):14, 2021. doi:10.5334/jors.309.\n",
+      "[5] Robert Timms, Scott G Marquis, Valentin Sulzer, Colin P. Please, and S Jonathan Chapman. Asymptotic Reduction of a Lithium-ion Pouch Cell Model. SIAM Journal on Applied Mathematics, 81(3):765–788, 2021. doi:10.1137/20M1336898.\n",
       "\n"
      ]
     }
@@ -165,7 +173,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.11.1"
+   "version": "3.9.16"
   },
   "vscode": {
    "interpreter": {

pybamm/models/submodels/electrolyte_diffusion/full_diffusion.py

@@ -89,17 +89,12 @@ def get_coupled_variables(self, variables):
     def set_rhs(self, variables):
         eps_c_e = variables["Porosity times concentration [mol.m-3]"]
         c_e = variables["Electrolyte concentration [mol.m-3]"]
-        T = variables["Cell temperature [K]"]
-        N_e_diffusion = variables["Electrolyte diffusion flux [mol.m-2.s-1]"]
-        N_e_convection = variables["Electrolyte convection flux [mol.m-2.s-1]"]
+        N_e = variables["Electrolyte flux [mol.m-2.s-1]"]
         div_Vbox = variables["Transverse volume-averaged acceleration [m.s-2]"]
 
         sum_s_a_j = variables["Sum of electrolyte reaction source terms [A.m-3]"]
-        sum_a_j = variables["Sum of volumetric interfacial current densities [A.m-3]"]
         sum_s_a_j.print_name = "aj"
-        source_terms = (sum_s_a_j - self.param.t_plus(c_e, T) * sum_a_j) / self.param.F
-
-        N_e = N_e_diffusion + N_e_convection
+        source_terms = sum_s_a_j / self.param.F
 
         self.rhs = {eps_c_e: -pybamm.div(N_e) + source_terms - c_e * div_Vbox}
 

tests/integration/test_models/standard_output_comparison.py

@@ -128,13 +128,13 @@ def test_all(self):
         self.compare("X-averaged positive electrode open-circuit potential [V]")
         self.compare("Voltage [V]")
         self.compare("X-averaged solid phase ohmic losses [V]")
-        self.compare("Negative electrode reaction overpotential [V]")
+        self.compare("Negative electrode reaction overpotential [V]", atol=1e-4)
         self.compare("Positive electrode reaction overpotential [V]")
         self.compare("Negative electrode potential [V]", atol=1e-5)
         self.compare("Positive electrode potential [V]")
         self.compare("Electrolyte potential [V]")
         # Currents
-        self.compare("Exchange current density [A.m-2]")
+        self.compare("Exchange current density [A.m-2]", atol=2e-3)
         self.compare("Negative electrode current density [A.m-2]", atol=1e-10)
         self.compare("Positive electrode current density [A.m-2]", atol=1e-10)
 

tests/integration/test_models/standard_output_tests.py

@@ -545,7 +545,10 @@ def test_conservation(self):
         diff = (
             self.c_e_tot(self.solution.t[1:]) - self.c_e_tot(self.solution.t[:-1])
         ) / self.c_e_tot(self.solution.t[:-1])
-        if self.model.options["surface form"] == "differential":
+        if self.model.options["surface form"] == "differential" or (
+            isinstance(self.model, pybamm.lithium_ion.DFN)
+            and self.model.options["surface form"] == "algebraic"
+        ):
             np.testing.assert_allclose(0, diff, atol=1e-4, rtol=1e-6)
         else:
             np.testing.assert_allclose(0, diff, atol=1e-14, rtol=1e-14)