Volume fractions and porosity that vary in y and z

CHANGELOG.md

@@ -9,6 +9,7 @@
 ## Features
 - Allow setting initial conditions from `y_slices` of a `Solution` object. ([#5257](https://github.com/pybamm-team/PyBaMM/pull/5257))
 - Added docstring to `FuzzyDict.copy` explaining its return value and behavior. ([#5242](https://github.com/pybamm-team/PyBaMM/pull/5242))
+- Porosity and active material fractions are now `FunctionParameters` of y and z, as well as x ([#5214](https://github.com/pybamm-team/PyBaMM/pull/5214))
 
 ## Bug fixes
 

docs/source/examples/notebooks/models/graded-electrodes.ipynb

@@ -22,13 +22,6 @@
      "text": [
       "Note: you may need to restart the kernel to use updated packages.\n"
      ]
-    },
-    {
-     "name": "stderr",
-     "output_type": "stream",
-     "text": [
-      "An NVIDIA GPU may be present on this machine, but a CUDA-enabled jaxlib is not installed. Falling back to cpu.\n"
-     ]
     }
    ],
    "source": [
@@ -59,7 +52,7 @@
    "source": [
     "We will vary the porosity in both electrodes and we will try three different scenarios: constant porosity, one where lower porosity occurs near the separator and one where lower porosity occurs near the current collector. All other parameters are kept constant. The varying porosity is defined to be linear centered around the default value and with a variation of $\\pm$ 10%.\n",
     "\n",
-    "We define the varying porosities and store them in a list so we can loop over when solving the model."
+    "We define the varying porosities and store them in a list so we can loop over when solving the model. Note that the porosities must be specified as functions of x, y, and z, even if y and z are not present in the model and not referred to at any other point, as in this example."
    ]
   },
   {
@@ -77,14 +70,14 @@
     "\n",
     "eps_ns = [\n",
     "    eps_n_0,\n",
-    "    lambda x: eps_n_0 * (1.1 - 0.2 * (x / L_n)),\n",
-    "    lambda x: eps_n_0 * (0.9 + 0.2 * (x / L_n)),\n",
+    "    lambda x, y, z: eps_n_0 * (1.1 - 0.2 * (x / L_n)),\n",
+    "    lambda x, y, z: eps_n_0 * (0.9 + 0.2 * (x / L_n)),\n",
     "]\n",
     "\n",
     "eps_ps = [\n",
     "    eps_p_0,\n",
-    "    lambda x: eps_p_0 * (0.9 - 0.2 / L_p * (L_n + L_s) + 0.2 * (x / L_p)),\n",
-    "    lambda x: eps_p_0 * (1.1 + 0.2 / L_p * (L_n + L_s) - 0.2 * (x / L_p)),\n",
+    "    lambda x, y, z: eps_p_0 * (0.9 - 0.2 / L_p * (L_n + L_s) + 0.2 * (x / L_p)),\n",
+    "    lambda x, y, z: eps_p_0 * (1.1 + 0.2 / L_p * (L_n + L_s) - 0.2 * (x / L_p)),\n",
     "]"
    ]
   },
@@ -132,12 +125,12 @@
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "99f847ca09da40cba550dd02dd8281a2",
+       "model_id": "9198b6c974714d07b0a68f1cbc42c988",
        "version_major": 2,
        "version_minor": 0
       },
       "text/plain": [
-       "interactive(children=(FloatSlider(value=0.0, description='t', max=673.9136958613059, step=6.7391369586130585),…"
+       "interactive(children=(FloatSlider(value=0.0, description='t', max=673.9107055793413, step=6.739107055793413), …"
       ]
      },
      "metadata": {},
@@ -146,7 +139,7 @@
     {
      "data": {
       "text/plain": [
-       "<pybamm.plotting.quick_plot.QuickPlot at 0x7fb504fb6f90>"
+       "<pybamm.plotting.quick_plot.QuickPlot at 0x77a5487643d0>"
       ]
      },
      "execution_count": 5,
@@ -176,18 +169,18 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 8,
+   "execution_count": 6,
    "metadata": {},
    "outputs": [
     {
      "data": {
       "application/vnd.jupyter.widget-view+json": {
-       "model_id": "829b68c6b3e04e0ebe5537daabec2278",
+       "model_id": "666bc13bf8cb4ee9b1a110cc828f025b",
        "version_major": 2,
        "version_minor": 0
       },
       "text/plain": [
-       "interactive(children=(FloatSlider(value=0.0, description='t', max=673.9136958613059, step=6.7391369586130585),…"
+       "interactive(children=(FloatSlider(value=0.0, description='t', max=673.9107055793413, step=6.739107055793413), …"
       ]
      },
      "metadata": {},
@@ -196,10 +189,10 @@
     {
      "data": {
       "text/plain": [
-       "<pybamm.plotting.quick_plot.QuickPlot at 0x7fb4fc650410>"
+       "<pybamm.plotting.quick_plot.QuickPlot at 0x77a5487674c0>"
       ]
      },
-     "execution_count": 8,
+     "execution_count": 6,
      "metadata": {},
      "output_type": "execute_result"
     }
@@ -218,7 +211,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 9,
+   "execution_count": 7,
    "metadata": {},
    "outputs": [
     {
@@ -232,10 +225,7 @@
       "[5] 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",
       "[6] Alan C. Hindmarsh. The PVODE and IDA algorithms. Technical Report, Lawrence Livermore National Lab., CA (US), 2000. doi:10.2172/802599.\n",
       "[7] Alan C. Hindmarsh, Peter N. Brown, Keith E. Grant, Steven L. Lee, Radu Serban, Dan E. Shumaker, and Carol S. Woodward. SUNDIALS: Suite of nonlinear and differential/algebraic equation solvers. ACM Transactions on Mathematical Software (TOMS), 31(3):363–396, 2005. doi:10.1145/1089014.1089020.\n",
-      "[8] Peyman Mohtat, Suhak Lee, Jason B Siegel, and Anna G Stefanopoulou. Towards better estimability of electrode-specific state of health: decoding the cell expansion. Journal of Power Sources, 427:101–111, 2019.\n",
-      "[9] 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",
-      "[10] Pauli Virtanen, Ralf Gommers, Travis E. Oliphant, Matt Haberland, Tyler Reddy, David Cournapeau, Evgeni Burovski, Pearu Peterson, Warren Weckesser, Jonathan Bright, and others. SciPy 1.0: fundamental algorithms for scientific computing in Python. Nature Methods, 17(3):261–272, 2020. doi:10.1038/s41592-019-0686-2.\n",
-      "[11] Andrew Weng, Jason B Siegel, and Anna Stefanopoulou. Differential voltage analysis for battery manufacturing process control. arXiv preprint arXiv:2303.07088, 2023.\n",
+      "[8] 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",
       "\n"
      ]
     }
@@ -254,7 +244,7 @@
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "venv",
+   "display_name": "Python 3 (ipykernel)",
    "language": "python",
    "name": "python3"
   },
@@ -268,9 +258,9 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.11.6"
+   "version": "3.10.12"
   }
  },
  "nbformat": 4,
- "nbformat_minor": 2
+ "nbformat_minor": 4
 }

examples/scripts/3d_examples/dried_out_pouch.py

@@ -0,0 +1,48 @@
+import numpy as np
+
+import pybamm
+
+model = pybamm.lithium_ion.DFN(
+    {"current collector": "potential pair", "dimensionality": 2}
+)
+param = pybamm.ParameterValues("Ecker2015")
+Ly = param["Electrode width [m]"]
+Lz = param["Electrode height [m]"]
+
+
+def _sigmoid(arg):
+    return (1 + np.tanh(arg)) / 2
+
+
+def _top_hat(arg, a, b, k=500):
+    return _sigmoid(k * (arg - a)) * _sigmoid(k * (b - arg))
+
+
+# Simulate drying out of the negative electrode edges by reducing porosity
+def eps_n(x, y, z):
+    return 0.329 * (
+        _top_hat(arg=y, a=Ly * 0.02, b=Ly * 0.98)
+        * _top_hat(arg=z, a=Lz * 0.02, b=Lz * 0.98)
+    )
+
+
+param.update({"Negative electrode porosity": eps_n})
+var_pts = {"x_n": 8, "x_s": 8, "x_p": 8, "r_n": 8, "r_p": 8, "y": 24, "z": 24}
+exp = pybamm.Experiment(
+    [
+        "Discharge at 1C until 2.7 V",
+        "Charge at 1C until 4.2 V",
+        "Hold at 4.2 V until C/20",
+    ]
+)
+sim = pybamm.Simulation(model, var_pts=var_pts, parameter_values=param, experiment=exp)
+sol = sim.solve()
+output_variables = [
+    "X-averaged negative electrode porosity",
+    "X-averaged negative particle surface stoichiometry",
+    "X-averaged negative electrode surface potential difference [V]",
+    "Current collector current density [A.m-2]",
+    "Current [A]",
+    "Voltage [V]",
+]
+plot = sol.plot(output_variables, variable_limits="tight", shading="auto")

src/pybamm/parameters/lithium_ion_parameters.py

@@ -233,8 +233,15 @@ def _set_parameters(self):
 
         if domain == "separator":
             x = pybamm.standard_spatial_vars.x_s
+            y = pybamm.PrimaryBroadcast(pybamm.standard_spatial_vars.y, "separator")
+            z = pybamm.PrimaryBroadcast(pybamm.standard_spatial_vars.z, "separator")
             self.epsilon_init = pybamm.FunctionParameter(
-                "Separator porosity", {"Through-cell distance (x) [m]": x}
+                "Separator porosity",
+                {
+                    "Through-cell distance (x) [m]": x,
+                    "Horizontal distance (y) [m]": y,
+                    "Vertical distance (z) [m]": z,
+                },
             )
             self.epsilon_inactive = 1 - self.epsilon_init
             return
@@ -248,6 +255,12 @@ def _set_parameters(self):
             auxiliary_domains={"secondary": "current collector"},
             coord_sys="cartesian",
         )
+        y = pybamm.PrimaryBroadcast(
+            pybamm.standard_spatial_vars.y, f"{domain} electrode"
+        )
+        z = pybamm.PrimaryBroadcast(
+            pybamm.standard_spatial_vars.z, f"{domain} electrode"
+        )
 
         # Macroscale geometry
         self.L_cc = self.geo.L_cc
@@ -267,7 +280,12 @@ def _set_parameters(self):
         )
         if main.options.electrode_types[domain] == "porous":
             self.epsilon_init = pybamm.FunctionParameter(
-                f"{Domain} electrode porosity", {"Through-cell distance (x) [m]": x}
+                f"{Domain} electrode porosity",
+                {
+                    "Through-cell distance (x) [m]": x,
+                    "Horizontal distance (y) [m]": y,
+                    "Vertical distance (z) [m]": z,
+                },
             )
             epsilon_s_tot = sum(phase.epsilon_s for phase in self.phase_params.values())
             self.epsilon_inactive = 1 - self.epsilon_init - epsilon_s_tot
@@ -414,6 +432,12 @@ def _set_parameters(self):
             auxiliary_domains={"secondary": "current collector"},
             coord_sys="cartesian",
         )
+        y = pybamm.PrimaryBroadcast(
+            pybamm.standard_spatial_vars.y, f"{domain} electrode"
+        )
+        z = pybamm.PrimaryBroadcast(
+            pybamm.standard_spatial_vars.z, f"{domain} electrode"
+        )
         r = pybamm.SpatialVariable(
             f"r_{domain[0]}",
             domain=[f"{domain} {self.phase_name}particle"],
@@ -438,7 +462,11 @@ def _set_parameters(self):
         # Particle properties
         self.epsilon_s = pybamm.FunctionParameter(
             f"{pref}{Domain} electrode active material volume fraction",
-            {"Through-cell distance (x) [m]": x},
+            {
+                "Through-cell distance (x) [m]": x,
+                "Horizontal distance (y) [m]": y,
+                "Vertical distance (z) [m]": z,
+            },
         )
         self.epsilon_s_av = pybamm.xyz_average(self.epsilon_s)
         self.c_max = pybamm.Parameter(