diff --git a/docs/_static/model_hierarchy.png b/docs/_static/model_hierarchy.png
new file mode 100644
index 0000000000..34411d183b
Binary files /dev/null and b/docs/_static/model_hierarchy.png differ
diff --git a/docs/source/api/models/submodels/electrolyte_conductivity/index.rst b/docs/source/api/models/submodels/electrolyte_conductivity/index.rst
index a723c32aab..8d0be2223a 100644
--- a/docs/source/api/models/submodels/electrolyte_conductivity/index.rst
+++ b/docs/source/api/models/submodels/electrolyte_conductivity/index.rst
@@ -7,6 +7,5 @@ Electrolyte Conductivity
     base_electrolyte_conductivity
     leading_order_conductivity
     composite_conductivity
-    integrated_conductivity
     full_conductivity
     surface_form/index
diff --git a/docs/source/api/models/submodels/electrolyte_conductivity/integrated_conductivity.rst b/docs/source/api/models/submodels/electrolyte_conductivity/integrated_conductivity.rst
deleted file mode 100644
index 144e1a19c7..0000000000
--- a/docs/source/api/models/submodels/electrolyte_conductivity/integrated_conductivity.rst
+++ /dev/null
@@ -1,7 +0,0 @@
-Integrated Model
-================
-
-.. autoclass:: pybamm.electrolyte_conductivity.Integrated
-    :members:
-
-.. footbibliography::
diff --git a/docs/source/user_guide/fundamentals/battery_models.md b/docs/source/user_guide/fundamentals/battery_models.md
index b2faff98dd..c552124d15 100644
--- a/docs/source/user_guide/fundamentals/battery_models.md
+++ b/docs/source/user_guide/fundamentals/battery_models.md
@@ -4,26 +4,190 @@ References for the battery models used in PyBaMM simulations can be found callin
 ```python
 pybamm.print_citations()
 ```
-However, a few papers are provided in this section for anyone interested in reading the theory
-behind the models before doing the tutorials.
+However, each model can be defined in slightly different ways depending on the paper. For this reason, in this page we state some of the standard battery models with the exact definition used in PyBaMM. The goal is not to provide a thorough explanation of the models, but just to state them. For more details, we strongly advise the reader to follow the references provided in this page.
 
-## Review Articles
+## Lithium-ion battery models
 
-[Review of physics-based lithium-ion battery models](https://doi.org/10.1088/2516-1083/ac7d31)
+The standard models for lithium-ion batteries can be classified in a hierarchical structure, from simple to complex. Here, as shown in the figure, we focus on the SPM, SPMe and DFN models. This hierarchy is discussed in a lot more detail in the review article by [Brosa Planella et al. (2022)](https://doi.org/10.1088/2516-1083/ac7d31), and also explains in detail how these models are derived.
 
-[Review of parameterisation and a novel database for Li-ion battery models](https://doi.org/10.1088/2516-1083/ac692c)
+<img title="Lithium-ion battery model hierarchy" alt="A diagram showing the hierarchy of lithium-ion battery models, from simple (SPM) to complex (DFN). The diagram also includes what each model captures and misses compared to the others." src="../../../_static/model_hierarchy.png">
 
-## Model References
+These models impose conservation of lithium and charge, so the variables we are interested in are the lithium concentration in the solid and electrolyte phases, $c_k$ and $c_\mathrm{e}$ respectively, and the electric potential $\phi_k$ and $\phi_\mathrm{e}$ in the solid and electrolyte phases, respectively. The subscript $k \in \{\mathrm{n}, \mathrm{p}\}$ denotes the negative and positive electrodes, respectively.
 
-### Lithium-Ion Batteries
+**Important remark:** these models account for the electrochemistry of the battery, and thus not include the thermal and degradation effects. These effects can be added to the electrochemical models as additional equations (and through the PyBaMM model options).
 
-[Doyle-Fuller-Newman model](https://doi.org/10.1149/1.2221597)
+The parameters involved in these models are defined in the following table (where parameters with $k \in \{\mathrm{n}, \mathrm{s}, \mathrm{p}\}$ may have different values for each electrode and separator):
 
-[Single particle model](https://doi.org/10.1149/2.0341915jes)
+| Symbol | Description | Units |
+|:------:|:-----------:|:-----:|
+| $R_k$ | Particle radius | m |
+| $L_k$ | Electrode/separator thickness | m |
+| $D_k (c_k)$ | Electrode particle diffusivity | m $^2$ s $^{-1}$ |
+| $\sigma_k$ | Electrode conductivity | S m $^{-1}$ |
+| $U_k (c_k)$ | Open-circuit potential | V |
+| $j_{0 k} (c_k, c_\mathrm{e})$ | Exchange current density | A m $^{-2}$ |
+| $\varepsilon_{\mathrm{s}, k}$ | Solid phase volume fraction | - |
+| $\varepsilon$ | Porosity | - |
+| $\mathcal{B}$ | Transport efficiency | - |
+| $D_\mathrm{e}(c_\mathrm{e})$ | Electrolyte diffusivity | m $^2$ s $^{-1}$ |
+| $\sigma_\mathrm{e}(c_\mathrm{e})$ | Electrolyte conductivity | S m $^{-1}$ |
+| $t^+(c_\mathrm{e})$ | Transference number | - |
+| $f_\mathrm{therm}(c_\mathrm{e})$ | Thermodynamic factor | - |
+| $c_\mathrm{e0}$ | Initial electrolyte concentration | mol m $^{-3}$ |
+| $T$ | Temperature | K |
+| $i_\mathrm{app}$ | Applied current density | A m $^{-2}$ |
 
+The surface area per unite volume $a_k$ is defined as $a_k = 3 \varepsilon_{\mathrm{s}, k}/ R_k$. The two physical constants appearing in the models are the Faraday constant is $F$, the gas constant $R$.
 
-### Lead-Acid Batteries
+Note that not all models use all of these parameters. For example, the SPM only uses a subset.
 
-[Isothermal porous-electrode model](https://doi.org/10.1149/2.0301910jes)
+### Single particle model (SPM)
+The SPM is the simplest of the three models considered here. It assumes that the active material particles behave similarly so they can be described by an averaged particle. The model considers two particles: one for the positive electrode and one for the negative electrode, where we model the lithium concentration as a function of space and time, denoted as $c_\mathrm{p}$ and $c_\mathrm{n}$ respectively. The model assumes spherical symmetry and the domains are $0 \leq r \leq R_k$. Two (independent) diffusion equations need to be solved for the two particles, and any additional quantities (such as the cell voltage) are computed from explicit expressions.
 
-[Leading-Order Quasi-Static model](https://doi.org/10.1149/2.0441908jes)
+
+In PyBaMM, the Single Particle Model is defined as follows:
+$$\begin{align}
+\frac{\partial c_\mathrm{n}}{\partial t} &= \frac{1}{r^2} \frac{\partial}{\partial r} \left (r^2 D_\mathrm{n} (c_\mathrm{n}) \frac{\partial c_\mathrm{n}}{\partial r} \right), & \text{ in } \quad 0 < r < R_\mathrm{n},\\
+\frac{\partial c_\mathrm{p}}{\partial t} &= \frac{1}{r^2} \frac{\partial}{\partial r} \left (r^2 D_\mathrm{p} (c_\mathrm{p}) \frac{\partial c_\mathrm{p}}{\partial r} \right), & \text{ in } \quad 0 < r < R_\mathrm{p},
+\end{align}$$
+
+with boundary conditions
+$$\begin{align}
+\frac{\partial c_\mathrm{n}}{\partial r} &= 0, & \text{ at } \quad r = 0,\\
+-D_\mathrm{p} (c_\mathrm{n}) \frac{\partial c_\mathrm{n}}{\partial r} &= \frac{i_\mathrm{app}(t)}{a_\mathrm{n} L_\mathrm{n} F}, & \text{ at } \quad r = R_\mathrm{n},\\
+\frac{\partial c_\mathrm{p}}{\partial r} &= 0, & \text{ at } \quad r = 0,\\
+-D_\mathrm{p} (c_\mathrm{p}) \frac{\partial c_\mathrm{p}}{\partial r} &= -\frac{i_\mathrm{app}(t)}{a_\mathrm{p} L_\mathrm{p} F}, & \text{ at } \quad r = R_\mathrm{p},
+\end{align}$$
+
+and initial conditions
+$$\begin{align}
+c_\mathrm{n}(r, 0) &= c_\mathrm{n0}(r), & \text{ at } \quad t = 0,\\
+c_\mathrm{p}(r, 0) &= c_\mathrm{p0}(r), & \text{ at } \quad t = 0.
+\end{align}
+$$
+
+The voltage can then be computed from $c_\mathrm{n}$ and $c_\mathrm{p}$ as
+$$\begin{align}
+V(t) &= U_\mathrm{p}(c_\mathrm{p}(R_\mathrm{p}, t)) - U_\mathrm{n}(c_\mathrm{n}(R_\mathrm{n}, t)) \\
+& \quad - \frac{2 R T}{F} \mathrm{arcsinh} \left(\frac{i_\mathrm{app}(t)}{a_\mathrm{n} L_\mathrm{n} j_{0\mathrm{n}}} \right) - \frac{2 R T}{F} \mathrm{arcsinh} \left(\frac{i_\mathrm{app}(t)}{a_\mathrm{p} L_\mathrm{p} j_{0\mathrm{p}}} \right),\nonumber
+\end{align}$$
+
+where $j_{0\mathrm{n}}(c_\mathrm{n})$ and $j_{0\mathrm{p}}(c_\mathrm{p})$ are the interfacial current densities at the negative and positive electrodes, respectively, and can be defined in the parameter set by the user.
+
+### Single particle model with electrolyte (SPMe)
+The SPMe enhances the SPM by including an electrolyte phase. This model is more complex than the SPM, but still relatively simple. In addition to the two particles, the model considers the electrolyte concentration $c_\mathrm{e}$ over the electrodes (and separator) thickness $0 \leq x \leq L$ (where $L = L_\mathrm{n} + L_\mathrm{s} + L_\mathrm{p}$).
+
+The equation for the electrolyte concentration is given by
+$$\begin{equation}
+\varepsilon\frac{\partial c_\mathrm{e}}{\partial t} = \frac{\partial}{\partial x} \left(D_\mathrm{e}(c_\mathrm{e}) \mathcal{B}(x) \frac{\partial c_\mathrm{e}}{\partial x} + \frac{1 - t^+(c_\mathrm{e})}{F} i_\mathrm{e}\right),
+\end{equation}$$
+with boundary conditions
+$$\begin{align}
+\frac{\partial c_\mathrm{e}}{\partial x} &= 0, & \text{ at } \quad x = 0 \quad \text{and} \quad x = L,\\
+\end{align}$$
+and initial condition
+$$\begin{align}
+c_\mathrm{e}(x, 0) &= c_\mathrm{e0}, & \text{ at } \quad t = 0.
+\end{align}$$
+
+The current in the electrolyte can be explicitly computed as
+$$\begin{equation}
+i_\mathrm{e}(x, t) = \begin{cases}
+\frac{i_\mathrm{app}(t)}{L_\mathrm{n}} x, & \text{ in } \quad 0 \leq x < L_\mathrm{n},\\
+i_\mathrm{app}(t), & \text{ in } \quad L_\mathrm{n} \leq x < L - L_\mathrm{p},\\
+\frac{i_\mathrm{app}(t)}{L_\mathrm{p}} (L - x), & \text{ in } \quad L - L_\mathrm{p} \leq x \leq L.
+\end{cases}
+\end{equation}$$
+
+From the electrolyte concentration, we can compute additional corrections to the SPM voltage expression. Then, the voltage can be computed as
+$$\begin{align}
+V(t) &= U_\mathrm{eq} + \eta_\mathrm{r} + \eta_\mathrm{e} + \Delta \phi_\mathrm{e} - \Delta \phi_\mathrm{s},
+\end{align}$$
+where each component is defined as
+$$\begin{align}
+U_\mathrm{eq} &= U_\mathrm{p}(c_\mathrm{p}(R_\mathrm{p}, t)) - U_\mathrm{n}(c_\mathrm{n}(R_\mathrm{n}, t)),\\
+\eta_\mathrm{r} &= -\frac{2 R T}{F} \left( \frac{1}{L_\mathrm{p}}\int_{L - L_\mathrm{p}}^L\mathrm{arcsinh} \left(\frac{i_\mathrm{app}(t)}{a_\mathrm{p} L_\mathrm{p} j_{0\mathrm{p}}} \right) \mathrm{d} x + \frac{1}{L_\mathrm{n}} \int_0^{L_\mathrm{n}}\mathrm{arcsinh} \left(\frac{i_\mathrm{app}(t)}{a_\mathrm{n} L_\mathrm{n} j_{0\mathrm{n}}} \right) \mathrm{d} x \right),\\
+\eta_\mathrm{e} &= \frac{2 R T}{F} \left( \frac{1}{L_\mathrm{p}}\int_{L - L_\mathrm{p}}^L \int_0^x (1 - t^+)f_\mathrm{therm} \frac{\partial \log c_\mathrm{e}(s, t)}{\partial s} \mathrm{d} s \mathrm{d} x + \frac{1}{L_\mathrm{n}} \int_0^{L_\mathrm{n}}\int_0^x (1 - t^+)f_\mathrm{therm} \frac{\partial \log c_\mathrm{e}(s, t)}{\partial s} \mathrm{d} s \mathrm{d} x \right),\\
+\Delta \phi_\mathrm{e} &= - \frac{1}{L_\mathrm{p}} \int_{L - L_\mathrm{p}}^L \int_0^x\frac{i_\mathrm{app}(t)}{\sigma_\mathrm{e}(c_\mathrm{e}(s, t)) \mathcal{B}(s)} \mathrm{d} s \mathrm{d} x + \frac{1}{L_\mathrm{n}} \int_0^{L_\mathrm{n}} \int_0^x\frac{i_\mathrm{app}(t)}{\sigma_\mathrm{e}(c_\mathrm{e}(s, t)) \mathcal{B}(s)} \mathrm{d} s \mathrm{d} x,\\
+\Delta \phi_\mathrm{s} &= -\frac{i_\mathrm{app}(t)}{3} \left( \frac{L_\mathrm{p}}{\sigma_\mathrm{p}} + \frac{L_\mathrm{n}}{\sigma_\mathrm{n}} \right),
+\end{align}$$
+where now $j_{0k}(c_k, c_\mathrm{e})$ depends both on the particle surface concentration and the electrolyte concentration.
+
+### Doyle-Fuller-Newman model (DFN)
+Finally, we present the DFN model, which is the most complex of the models considered here. The model still considers the multiscale domain of the particles $0 \leq r \leq R_k$ and the electrolyte $0 \leq x \leq L$, but now the model allows for a different particle at each point $x$ of the macroscale. Then, the model requires solving for the concentrations and potentials in both electrodes and the electrolyte: $c_k(r, x, t)$, $c_\mathrm{e}(x, t)$, $\phi_k(x, t)$, $\phi_\mathrm{e}(x, t)$.
+
+The concentration in the particles is described at the microscale by
+$$\begin{align}
+\frac{\partial c_k}{\partial t} &= \frac{1}{r^2} \frac{\partial}{\partial r} \left (r^2 D_k (c_k) \frac{\partial c_k}{\partial r} \right), & \text{ in } \quad 0 < r < R_k,\\
+\end{align}$$
+with boundary and initial conditions
+$$\begin{align}
+\frac{\partial c_k}{\partial r} &= 0, & \text{ at } \quad r = 0,\\
+-D_k (c_k) \frac{\partial c_k}{\partial r} &= \frac{J_k}{a_k F}, & \text{ at } \quad r = R_k,\\
+c_k(r, x, 0) &= c_{k0}(r, x), & \text{ at } \quad t = 0.
+\end{align}$$
+
+The concentration in the electrolyte is described by
+$$\begin{equation}
+\varepsilon\frac{\partial c_\mathrm{e}}{\partial t} = \frac{\partial}{\partial x} \left(D_\mathrm{e}(c_\mathrm{e}) \mathcal{B}(x) \frac{\partial c_\mathrm{e}}{\partial x} + \frac{1 - t^+(c_\mathrm{e})}{F} i_\mathrm{e}\right), \quad \text{ in } \quad 0 < x < L,
+\end{equation}$$
+with boundary and initial conditions
+$$\begin{align}
+\frac{\partial c_\mathrm{e}}{\partial x} &= 0, & \text{ at } \quad x = 0 \quad \text{and} \quad x = L,\\
+c_\mathrm{e}(x, 0) &= c_\mathrm{e0}, & \text{ at } \quad t = 0.
+\end{align}$$
+
+The potential in the electrodes is described by
+$$\begin{align}
+\frac{\partial i_k}{\partial x} = -J_k,
+\end{align}$$
+with
+$$\begin{align}
+i_k = -\sigma_k \frac{\partial \phi_k}{\partial x},
+\end{align}$$
+which for $k = \mathrm{n}$ is defined in $x \in [0, L_\mathrm{n}]$, and for $k = \mathrm{p}$ is defined in $x \in [L - L_\mathrm{p}, L]$. The boundary conditions are
+$$\begin{align}
+\phi_\mathrm{n} &= 0, & \text{ at } \quad x = 0,\\
+i_\mathrm{n} &= 0, & \text{ at } \quad x = L_\mathrm{n},\\
+i_\mathrm{p} &= 0, & \text{ at } \quad x = L - L_\mathrm{p},\\
+i_\mathrm{p} &= i_\mathrm{app}(t), & \text{ at } \quad x = L.
+\end{align}$$
+
+Finally, the potential in the electrolyte is described by
+$$\begin{align}
+\frac{\partial i_\mathrm{e}}{\partial x} = J, & \quad \text{ in } \quad 0 < x < L,
+\end{align}$$
+with
+$$\begin{align}
+i_\mathrm{e} = -\sigma_\mathrm{e} \mathcal{B}(x) \left(\frac{\partial \phi_\mathrm{e}}{\partial x} - 2 (1 - t^+) f_\mathrm{therm} \frac{R T}{F} \frac{\partial \log c_\mathrm{e}}{\partial x}\right),
+\end{align}$$
+and boundary conditions
+$$\begin{align}
+i_\mathrm{e} &= 0, & \text{ at } \quad x = 0 \quad \text{ and } \quad x = L.
+\end{align}$$
+
+The intercalation current density $J$ is defined as
+$$\begin{align}
+J = \begin{cases}
+J_\mathrm{n}, & \text{ in } \quad 0 \leq x \leq L_\mathrm{n},\\
+0, & \text{ in } \quad L_\mathrm{n} < x \leq L - L_\mathrm{p},\\
+J_\mathrm{p}, & \text{ in } \quad L - L_\mathrm{p} < x \leq L,
+\end{cases}
+\end{align}$$
+where $J_\mathrm{n}$ and $J_\mathrm{p}$ are the intercalation current densities at the negative and positive electrodes, respectively, defined as
+$$\begin{align}
+J_k = -a_k j_{0 k} \sinh \left( \frac{F}{2RT} \left(\phi_k - \phi_\mathrm{e} - U_k(c_k(R_k, t)) \right) \right).
+\end{align}$$
+
+## Further references
+
+1. M. Doyle, T.F. Fuller, J. Newman, [Modeling of Galvanostatic Charge and Discharge of the Lithium/Polymer/Insertion Cell](https://doi.org/10.1149/1.2221597), Journal of The Electrochemical Society 140 (1993) 1526–1533.
+
+
+2. T.F. Fuller, M. Doyle, J. Newman, [Simulation and Optimization of the Dual Lithium Ion Insertion Cell](https://doi.org/10.1149/1.2054684), Journal of The Electrochemical Society 141 (1994) 1–10.
+
+3. S.G. Marquis, V. Sulzer, R. Timms, C.P. Please, S.J. Chapman, [An Asymptotic Derivation of a Single Particle Model with Electrolyte](https://doi.org/10.1149/2.0341915jes), Journal of The Electrochemical Society 166 (2019) A3693–A3706.
+
+4. F. Brosa Planella, W. Ai, A.M. Boyce, A. Ghosh, I. Korotkin, S. Sahu, V. Sulzer, R. Timms, T.G. Tranter, M. Zyskin, S.J. Cooper, J.S. Edge, J.M. Foster, M. Marinescu, B. Wu, G. Richardson, [A continuum of physics-based lithium-ion battery models reviewed](https://doi.org/10.1088/2516-1083/ac7d31), Prog. Energy 4 (2022) 042003.
+
+5. F. Brosa Planella, W.D. Widanage, [A Single Particle model with electrolyte and side reactions for degradation of lithium-ion batteries](https://doi.org/10.1016/j.apm.2022.12.009), Applied Mathematical Modelling 121 (2023) 586–610.
diff --git a/src/pybamm/expression_tree/unary_operators.py b/src/pybamm/expression_tree/unary_operators.py
index aa90fd6f4c..5e898467bc 100644
--- a/src/pybamm/expression_tree/unary_operators.py
+++ b/src/pybamm/expression_tree/unary_operators.py
@@ -771,7 +771,7 @@ class IndefiniteIntegral(BaseIndefiniteIntegral):
     A node in the expression tree representing an indefinite integral operator.
 
     .. math::
-        I = \\int_{x_\text{min}}^{x}\\!f(u)\\,du
+        I = \\int_{x_{\text{min}}}^{x}\\!f(u)\\,du
 
     where :math:`u\\in\\text{domain}` which can represent either a
     spatial or temporal variable.
diff --git a/src/pybamm/models/full_battery_models/base_battery_model.py b/src/pybamm/models/full_battery_models/base_battery_model.py
index a445f47d19..4ad0b858c2 100644
--- a/src/pybamm/models/full_battery_models/base_battery_model.py
+++ b/src/pybamm/models/full_battery_models/base_battery_model.py
@@ -63,8 +63,7 @@ class BatteryModelOptions(pybamm.FuzzyDict):
                 Sets the dimension of the current collector problem. Can be 0
                 (default), 1 or 2.
             * "electrolyte conductivity" : str
-                Can be "default" (default), "full", "leading order", "composite" or
-                "integrated".
+                Can be "default" (default), "full", "leading order", "composite".
             * "exchange-current density" : str
                 Sets the model for the exchange-current density. Can be "single"
                 (default) or "current sigmoid". A 2-tuple can be provided for different
@@ -243,7 +242,6 @@ def __init__(self, extra_options):
                 "full",
                 "leading order",
                 "composite",
-                "integrated",
             ],
             "exchange-current density": ["single", "current sigmoid"],
             "heat of mixing": ["false", "true"],
@@ -943,7 +941,6 @@ def options(self, extra_options):
             if options["electrolyte conductivity"] not in [
                 "default",
                 "composite",
-                "integrated",
             ]:
                 raise pybamm.OptionError(
                     "electrolyte conductivity '{}' not suitable for SPMe".format(
diff --git a/src/pybamm/models/full_battery_models/lithium_ion/spme.py b/src/pybamm/models/full_battery_models/lithium_ion/spme.py
index 6e0784ff60..a07b8a82cf 100644
--- a/src/pybamm/models/full_battery_models/lithium_ion/spme.py
+++ b/src/pybamm/models/full_battery_models/lithium_ion/spme.py
@@ -61,12 +61,6 @@ def set_electrolyte_potential_submodel(self):
                         self.param, options=self.options
                     )
                 )
-            elif self.options["electrolyte conductivity"] == "integrated":
-                self.submodels["electrolyte conductivity"] = (
-                    pybamm.electrolyte_conductivity.Integrated(
-                        self.param, options=self.options
-                    )
-                )
         if self.options["surface form"] == "false":
             surf_model = surf_form.Explicit
         elif self.options["surface form"] == "differential":
diff --git a/src/pybamm/models/submodels/electrolyte_conductivity/__init__.py b/src/pybamm/models/submodels/electrolyte_conductivity/__init__.py
index 48b3232d44..686606e64d 100644
--- a/src/pybamm/models/submodels/electrolyte_conductivity/__init__.py
+++ b/src/pybamm/models/submodels/electrolyte_conductivity/__init__.py
@@ -2,7 +2,6 @@
 from .leading_order_conductivity import LeadingOrder
 from .composite_conductivity import Composite
 from .full_conductivity import Full
-from .integrated_conductivity import Integrated
 
 from . import surface_potential_form
 
diff --git a/src/pybamm/models/submodels/electrolyte_conductivity/composite_conductivity.py b/src/pybamm/models/submodels/electrolyte_conductivity/composite_conductivity.py
index d6c7ea6473..aa38bfa5ce 100644
--- a/src/pybamm/models/submodels/electrolyte_conductivity/composite_conductivity.py
+++ b/src/pybamm/models/submodels/electrolyte_conductivity/composite_conductivity.py
@@ -24,15 +24,13 @@ class Composite(BaseElectrolyteConductivity):
     def __init__(self, param, domain=None, options=None):
         super().__init__(param, domain, options=options)
 
-    def _higher_order_macinnes_function(self, x):
-        "Function to differentiate between composite and first-order models"
+    def _derivative_macinnes_function(self, x):
+        "Compute the derivative of the MacInnes function."
         tol = pybamm.settings.tolerances["macinnes__c_e"]
         x = pybamm.maximum(x, tol)
-        return pybamm.log(x)
+        return 1 / x
 
     def get_coupled_variables(self, variables):
-        c_e_av = variables["X-averaged electrolyte concentration [mol.m-3]"]
-
         i_boundary_cc = variables["Current collector current density [A.m-2]"]
         if self.options.electrode_types["negative"] == "porous":
             c_e_n = variables["Negative electrolyte concentration [mol.m-3]"]
@@ -40,134 +38,128 @@ def get_coupled_variables(self, variables):
                 "X-averaged negative electrode surface potential difference [V]"
             ]
             phi_s_n_av = variables["X-averaged negative electrode potential [V]"]
-            tor_n_av = variables["X-averaged negative electrolyte transport efficiency"]
+            tor_n = variables["Negative electrolyte transport efficiency"]
 
         c_e_s = variables["Separator electrolyte concentration [mol.m-3]"]
         c_e_p = variables["Positive electrolyte concentration [mol.m-3]"]
 
-        tor_s_av = variables["X-averaged separator electrolyte transport efficiency"]
-        tor_p_av = variables["X-averaged positive electrolyte transport efficiency"]
+        tor_s = variables["Separator electrolyte transport efficiency"]
+        tor_p = variables["Positive electrolyte transport efficiency"]
 
         T_av = variables["X-averaged cell temperature [K]"]
         T_av_s = pybamm.PrimaryBroadcast(T_av, "separator")
         T_av_p = pybamm.PrimaryBroadcast(T_av, "positive electrode")
 
-        RT_F_av = self.param.R * T_av / self.param.F
         RT_F_av_s = self.param.R * T_av_s / self.param.F
         RT_F_av_p = self.param.R * T_av_p / self.param.F
 
         L_n = self.param.n.L
-        L_s = self.param.s.L
         L_p = self.param.p.L
         L_x = self.param.L_x
         x_s = pybamm.standard_spatial_vars.x_s
         x_p = pybamm.standard_spatial_vars.x_p
-
-        # bulk conductivities
-        kappa_s_av = self.param.kappa_e(c_e_av, T_av) * tor_s_av
-        kappa_p_av = self.param.kappa_e(c_e_av, T_av) * tor_p_av
-
-        chi_av = self.param.chi(c_e_av, T_av)
-        chi_av_s = pybamm.PrimaryBroadcast(chi_av, "separator")
-        chi_av_p = pybamm.PrimaryBroadcast(chi_av, "positive electrode")
+        x_p_edge = pybamm.standard_spatial_vars.x_p_edge
 
         # electrolyte current
         if self.options.electrode_types["negative"] == "planar":
-            i_e_n = None
+            # i_e_n = None
+            i_e_n_edge = None
         else:
             x_n = pybamm.standard_spatial_vars.x_n
-            chi_av_n = pybamm.PrimaryBroadcast(chi_av, "negative electrode")
+            x_n_edge = pybamm.standard_spatial_vars.x_n_edge
             T_av_n = pybamm.PrimaryBroadcast(T_av, "negative electrode")
             RT_F_av_n = self.param.R * T_av_n / self.param.F
-            kappa_n_av = self.param.kappa_e(c_e_av, T_av) * tor_n_av
-            i_e_n = i_boundary_cc * x_n / L_n
-        i_e_s = pybamm.PrimaryBroadcast(i_boundary_cc, "separator")
-        i_e_p = i_boundary_cc * (L_x - x_p) / L_p
-        i_e = pybamm.concatenation(i_e_n, i_e_s, i_e_p)
-
+            # i_e_n = i_boundary_cc * x_n / L_n
+            i_e_n_edge = i_boundary_cc * x_n_edge / L_n
+        # i_e_s = pybamm.PrimaryBroadcast(i_boundary_cc, "separator")
+        i_e_s_edge = pybamm.PrimaryBroadcastToEdges(i_boundary_cc, "separator")
+        # i_e_p = i_boundary_cc * (L_x - x_p) / L_p
+        i_e_p_edge = i_boundary_cc * (L_x - x_p_edge) / L_p
+        i_e = pybamm.concatenation(i_e_n_edge, i_e_s_edge, i_e_p_edge)
         phi_e_dict = {}
 
         # electrolyte potential
         if self.options.electrode_types["negative"] == "planar":
-            phi_e_li = variables["Lithium metal interface electrolyte potential [V]"]
-            c_e_n = pybamm.boundary_value(c_e_s, "left")
-            phi_e_const = (
-                phi_e_li
-                - chi_av
-                * RT_F_av
-                * self._higher_order_macinnes_function(c_e_n / c_e_av)
-                + (i_boundary_cc / kappa_s_av) * L_n
-            )
+            phi_e_const = variables["Lithium metal interface electrolyte potential [V]"]
+            eta_c_n = pybamm.Scalar(0)
+            delta_phi_e_n = pybamm.Scalar(0)
         else:
+            # need to have current evaluated on edge as integral will make it
+            # evaluate on node
+
+            eta_c_n = -RT_F_av_n * pybamm.IndefiniteIntegral(
+                self.param.chi(c_e_n, T_av_n)
+                * self._derivative_macinnes_function(c_e_n)
+                * pybamm.grad(c_e_n),
+                x_n,
+            )
+
+            delta_phi_e_n = pybamm.IndefiniteIntegral(
+                i_e_n_edge / (self.param.kappa_e(c_e_n, T_av_n) * tor_n), x_n
+            )
             phi_e_const = (
                 -delta_phi_n_av
                 + phi_s_n_av
-                - (
-                    chi_av
-                    * RT_F_av
-                    * pybamm.x_average(
-                        self._higher_order_macinnes_function(c_e_n / c_e_av)
-                    )
-                )
-                - ((i_boundary_cc * L_n) * (1 / (3 * kappa_n_av) - 1 / kappa_s_av))
+                + pybamm.x_average(eta_c_n)
+                + pybamm.x_average(delta_phi_e_n)
             )
 
-            phi_e_n = (
-                phi_e_const
-                + (
-                    chi_av_n
-                    * RT_F_av_n
-                    * self._higher_order_macinnes_function(c_e_n / c_e_av)
-                )
-                - (i_boundary_cc / kappa_n_av) * (x_n**2 - L_n**2) / (2 * L_n)
-                - i_boundary_cc * L_n / kappa_s_av
-            )
+            phi_e_n = phi_e_const - (delta_phi_e_n + eta_c_n)
             phi_e_dict["negative electrode"] = phi_e_n
 
+        eta_c_s = -RT_F_av_s * pybamm.IndefiniteIntegral(
+            self.param.chi(c_e_s, T_av_s)
+            * self._derivative_macinnes_function(c_e_s)
+            * pybamm.grad(c_e_s),
+            x_s,
+        )
+        delta_phi_e_s = pybamm.IndefiniteIntegral(
+            i_e_s_edge / (self.param.kappa_e(c_e_s, T_av_s) * tor_s), x_s
+        )
+
         phi_e_s = (
             phi_e_const
-            + (
-                chi_av_s
-                * RT_F_av_s
-                * self._higher_order_macinnes_function(c_e_s / c_e_av)
-            )
-            - (i_boundary_cc / kappa_s_av) * x_s
+            - pybamm.boundary_value(delta_phi_e_n, "right")
+            - pybamm.boundary_value(eta_c_n, "right")
+            - (delta_phi_e_s + eta_c_s)
+        )
+
+        eta_c_p = -RT_F_av_p * pybamm.IndefiniteIntegral(
+            self.param.chi(c_e_p, T_av_p)
+            * self._derivative_macinnes_function(c_e_p)
+            * pybamm.grad(c_e_p),
+            x_p,
+        )
+        delta_phi_e_p = pybamm.IndefiniteIntegral(
+            i_e_p_edge / (self.param.kappa_e(c_e_p, T_av_p) * tor_p), x_p
         )
 
         phi_e_p = (
             phi_e_const
-            + (
-                chi_av_p
-                * RT_F_av_p
-                * self._higher_order_macinnes_function(c_e_p / c_e_av)
-            )
-            - (i_boundary_cc / kappa_p_av)
-            * (x_p * (2 * L_x - x_p) + L_p**2 - L_x**2)
-            / (2 * L_p)
-            - i_boundary_cc * (L_x - L_p) / kappa_s_av
+            - pybamm.boundary_value(delta_phi_e_n, "right")
+            - pybamm.boundary_value(eta_c_n, "right")
+            - pybamm.boundary_value(delta_phi_e_s, "right")
+            - pybamm.boundary_value(eta_c_s, "right")
+            - (delta_phi_e_p + eta_c_p)
         )
 
         phi_e_dict["separator"] = phi_e_s
         phi_e_dict["positive electrode"] = phi_e_p
 
         # concentration overpotential
-        macinnes_c_e_p = pybamm.x_average(
-            self._higher_order_macinnes_function(c_e_p / c_e_av)
+        eta_c_av = (
+            -pybamm.x_average(eta_c_p)
+            + pybamm.x_average(eta_c_n)
+            - pybamm.boundary_value(eta_c_s, "right")
+            - pybamm.boundary_value(eta_c_n, "right")
         )
-        if self.options.electrode_types["negative"] == "planar":
-            macinnes_c_e_n = 0
-            ohmic_n = 0
-        else:
-            macinnes_c_e_n = pybamm.x_average(
-                self._higher_order_macinnes_function(c_e_n / c_e_av)
-            )
-            ohmic_n = L_n / (3 * kappa_n_av)
-
-        eta_c_av = chi_av * RT_F_av * (macinnes_c_e_p - macinnes_c_e_n)
 
         # average electrolyte ohmic losses
-        delta_phi_e_av = -(i_boundary_cc) * (
-            ohmic_n + L_s / kappa_s_av + L_p / (3 * kappa_p_av)
+        delta_phi_e_av = (
+            -pybamm.x_average(delta_phi_e_p)
+            + pybamm.x_average(delta_phi_e_n)
+            - pybamm.boundary_value(delta_phi_e_s, "right")
+            - pybamm.boundary_value(delta_phi_e_n, "right")
         )
 
         variables.update(self._get_standard_potential_variables(phi_e_dict))
@@ -178,3 +170,42 @@ def get_coupled_variables(self, variables):
         i_e.print_name = "i_e"
 
         return variables
+
+    def set_boundary_conditions(self, variables):
+        # Define boundary conditions for electrolyte potential
+        super().set_boundary_conditions(variables)
+
+        # Define the boundary conditions for electrolyte concentration and potential
+        # so the gradients can be computed
+        if self.options.electrode_types["negative"] == "planar":
+            domains = ["Separator", "Positive"]
+        else:
+            domains = ["Negative", "Separator", "Positive"]
+
+        var_names = [
+            # " electrode porosity times concentration [mol.m-3]",
+            " electrolyte concentration [mol.m-3]",
+            # " electrolyte potential [V]",
+        ]
+
+        for name in var_names:
+            for domain in domains:
+                full_var_name = domain + name
+
+                # Fix inconcsistency in naming, where porosity times concentration is
+                # called "electrode porosity times concentration"
+                full_var_name = full_var_name.replace(
+                    "Separator electrode", "Separator"
+                )
+                var = variables[full_var_name]
+                self.boundary_conditions.update(
+                    {
+                        var: {
+                            "left": (pybamm.boundary_gradient(var, "left"), "Neumann"),
+                            "right": (
+                                pybamm.boundary_gradient(var, "right"),
+                                "Neumann",
+                            ),
+                        },
+                    }
+                )
diff --git a/src/pybamm/models/submodels/electrolyte_conductivity/integrated_conductivity.py b/src/pybamm/models/submodels/electrolyte_conductivity/integrated_conductivity.py
deleted file mode 100644
index 2250d99f6d..0000000000
--- a/src/pybamm/models/submodels/electrolyte_conductivity/integrated_conductivity.py
+++ /dev/null
@@ -1,164 +0,0 @@
-#
-# Composite electrolyte potential employing integrated Stefan-Maxwell
-#
-import pybamm
-from .base_electrolyte_conductivity import BaseElectrolyteConductivity
-
-
-class Integrated(BaseElectrolyteConductivity):
-    """
-    Integrated model for conservation of charge in the electrolyte derived from
-    integrating the Stefan-Maxwell constitutive equations, from
-    :footcite:t:`BrosaPlanella2021`.
-
-    Parameters
-    ----------
-    param : parameter class
-        The parameters to use for this submodel
-    domain : str, optional
-        The domain in which the model holds
-    options : dict, optional
-        A dictionary of options to be passed to the model.
-
-    """
-
-    def __init__(self, param, domain=None, options=None):
-        super().__init__(param, domain, options=options)
-        pybamm.citations.register("BrosaPlanella2021")
-
-    def _higher_order_macinnes_function(self, x):
-        tol = pybamm.settings.tolerances["macinnes__c_e"]
-        x = pybamm.maximum(x, tol)
-        return pybamm.log(x)
-
-    def get_coupled_variables(self, variables):
-        c_e_av = variables["X-averaged electrolyte concentration [mol.m-3]"]
-
-        i_boundary_cc = variables["Current collector current density [A.m-2]"]
-        c_e_n = variables["Negative electrolyte concentration [mol.m-3]"]
-        c_e_s = variables["Separator electrolyte concentration [mol.m-3]"]
-        c_e_p = variables["Positive electrolyte concentration [mol.m-3]"]
-        c_e_n0 = pybamm.boundary_value(c_e_n, "left")
-
-        delta_phi_n_av = variables[
-            "X-averaged negative electrode surface potential difference [V]"
-        ]
-        phi_s_n_av = variables["X-averaged negative electrode potential [V]"]
-
-        tor_n = variables["Negative electrolyte transport efficiency"]
-        tor_s = variables["Separator electrolyte transport efficiency"]
-        tor_p = variables["Positive electrolyte transport efficiency"]
-
-        T_av = variables["X-averaged cell temperature [K]"]
-        T_av_n = pybamm.PrimaryBroadcast(T_av, "negative electrode")
-        T_av_s = pybamm.PrimaryBroadcast(T_av, "separator")
-        T_av_p = pybamm.PrimaryBroadcast(T_av, "positive electrode")
-
-        RT_F_av = self.param.R * T_av / self.param.F
-        RT_F_av_n = self.param.R * T_av_n / self.param.F
-        RT_F_av_s = self.param.R * T_av_s / self.param.F
-        RT_F_av_p = self.param.R * T_av_p / self.param.F
-
-        L_n = self.param.n.L
-        L_p = self.param.p.L
-        L_x = self.param.L_x
-        x_n = pybamm.standard_spatial_vars.x_n
-        x_s = pybamm.standard_spatial_vars.x_s
-        x_p = pybamm.standard_spatial_vars.x_p
-        x_n_edge = pybamm.standard_spatial_vars.x_n_edge
-        x_p_edge = pybamm.standard_spatial_vars.x_p_edge
-
-        chi_av = self.param.chi(c_e_av, T_av)
-        chi_av_n = pybamm.PrimaryBroadcast(chi_av, "negative electrode")
-        chi_av_s = pybamm.PrimaryBroadcast(chi_av, "separator")
-        chi_av_p = pybamm.PrimaryBroadcast(chi_av, "positive electrode")
-
-        # electrolyte current
-        i_e_n = i_boundary_cc * x_n / L_n
-        i_e_s = pybamm.PrimaryBroadcast(i_boundary_cc, "separator")
-        i_e_p = i_boundary_cc * (L_x - x_p) / L_p
-        i_e = pybamm.concatenation(i_e_n, i_e_s, i_e_p)
-
-        i_e_n_edge = i_boundary_cc * x_n_edge / L_n
-        i_e_s_edge = pybamm.PrimaryBroadcastToEdges(i_boundary_cc, "separator")
-        i_e_p_edge = i_boundary_cc * (L_x - x_p_edge) / L_p
-
-        # electrolyte potential
-        indef_integral_n = pybamm.IndefiniteIntegral(
-            i_e_n_edge / (self.param.kappa_e(c_e_n, T_av_n) * tor_n), x_n
-        )
-        indef_integral_s = pybamm.IndefiniteIntegral(
-            i_e_s_edge / (self.param.kappa_e(c_e_s, T_av_s) * tor_s), x_s
-        )
-        indef_integral_p = pybamm.IndefiniteIntegral(
-            i_e_p_edge / (self.param.kappa_e(c_e_p, T_av_p) * tor_p), x_p
-        )
-
-        integral_n = indef_integral_n
-        integral_s = indef_integral_s + pybamm.boundary_value(integral_n, "right")
-        integral_p = indef_integral_p + pybamm.boundary_value(integral_s, "right")
-
-        phi_e_const = (
-            -delta_phi_n_av
-            + phi_s_n_av
-            - (
-                chi_av
-                * RT_F_av
-                * pybamm.x_average(self._higher_order_macinnes_function(c_e_n / c_e_n0))
-            )
-            + pybamm.x_average(integral_n)
-        )
-
-        phi_e_n = (
-            phi_e_const
-            + (
-                chi_av_n
-                * RT_F_av_n
-                * self._higher_order_macinnes_function(c_e_n / c_e_n0)
-            )
-            - integral_n
-        )
-
-        phi_e_s = (
-            phi_e_const
-            + (
-                chi_av_s
-                * RT_F_av_s
-                * self._higher_order_macinnes_function(c_e_s / c_e_n0)
-            )
-            - integral_s
-        )
-
-        phi_e_p = (
-            phi_e_const
-            + (
-                chi_av_p
-                * RT_F_av_p
-                * self._higher_order_macinnes_function(c_e_p / c_e_n0)
-            )
-            - integral_p
-        )
-
-        # concentration overpotential
-        eta_c_av = (
-            chi_av
-            * RT_F_av
-            * (
-                pybamm.x_average(self._higher_order_macinnes_function(c_e_p / c_e_av))
-                - pybamm.x_average(self._higher_order_macinnes_function(c_e_n / c_e_av))
-            )
-        )
-
-        # average electrolyte ohmic losses
-        delta_phi_e_av = -(pybamm.x_average(integral_p) - pybamm.x_average(integral_n))
-
-        phi_e_dict = {
-            "negative electrode": phi_e_n,
-            "separator": phi_e_s,
-            "positive electrode": phi_e_p,
-        }
-        variables.update(self._get_standard_potential_variables(phi_e_dict))
-        variables.update(self._get_standard_current_variables(i_e))
-        variables.update(self._get_split_overpotential(eta_c_av, delta_phi_e_av))
-
-        return variables
diff --git a/src/pybamm/models/submodels/thermal/base_thermal.py b/src/pybamm/models/submodels/thermal/base_thermal.py
index 42d90f1bcf..85eb076fda 100644
--- a/src/pybamm/models/submodels/thermal/base_thermal.py
+++ b/src/pybamm/models/submodels/thermal/base_thermal.py
@@ -102,8 +102,11 @@ def _get_standard_coupled_variables(self, variables):
                 i_e_n, i_e_s, i_e_p = i_e.orphans
                 phi_e_n = variables["Negative electrolyte potential [V]"]
                 Q_ohm_e_n = -pybamm.inner(i_e_n, pybamm.grad(phi_e_n))
+                # Q_ohm_e_n = -pybamm.inner(i_e_n, i_e_n)
             Q_ohm_e_s = -pybamm.inner(i_e_s, pybamm.grad(phi_e_s))
+            # Q_ohm_e_s = -pybamm.inner(i_e_s, i_e_s)
             Q_ohm_e_p = -pybamm.inner(i_e_p, pybamm.grad(phi_e_p))
+            # Q_ohm_e_p = -pybamm.inner(i_e_p, i_e_p)
             Q_ohm_e = pybamm.concatenation(Q_ohm_e_n, Q_ohm_e_s, Q_ohm_e_p)
         else:
             # else compute using i_e across all domains
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_spme.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_spme.py
index 6360b513c7..66345614e7 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_spme.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/test_spme.py
@@ -10,7 +10,3 @@ class TestSPMe(BaseIntegrationTestLithiumIon):
     @pytest.fixture(autouse=True)
     def setup(self):
         self.model = pybamm.lithium_ion.SPMe
-
-    def test_integrated_conductivity(self):
-        options = {"electrolyte conductivity": "integrated"}
-        self.run_basic_processing_test(options)
diff --git a/tests/unit/test_citations.py b/tests/unit/test_citations.py
index c87912490f..2c78a341a2 100644
--- a/tests/unit/test_citations.py
+++ b/tests/unit/test_citations.py
@@ -227,16 +227,6 @@ def test_subramanian_2005(self):
         assert "Subramanian2005" in citations._papers_to_cite
         assert "Subramanian2005" in citations._citation_tags.keys()
 
-    def test_brosaplanella_2021(self):
-        # Test that calling relevant bits of code adds the right paper to citations
-        citations = pybamm.citations
-
-        citations._reset()
-        assert "BrosaPlanella2021" not in citations._papers_to_cite
-        pybamm.electrolyte_conductivity.Integrated(None)
-        assert "BrosaPlanella2021" in citations._papers_to_cite
-        assert "BrosaPlanella2021" in citations._citation_tags.keys()
-
     def test_brosaplanella_2022(self):
         # Test that calling relevant bits of code adds the right paper to citations
         citations = pybamm.citations
diff --git a/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py b/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py
index 7dcfccdb66..4a7dbb1788 100644
--- a/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py
+++ b/tests/unit/test_models/test_full_battery_models/test_base_battery_model.py
@@ -24,7 +24,7 @@
 'current collector': 'uniform' (possible: ['uniform', 'potential pair', 'potential pair quite conductive'])
 'diffusivity': 'single' (possible: ['single', 'current sigmoid'])
 'dimensionality': 0 (possible: [0, 1, 2])
-'electrolyte conductivity': 'default' (possible: ['default', 'full', 'leading order', 'composite', 'integrated'])
+'electrolyte conductivity': 'default' (possible: ['default', 'full', 'leading order', 'composite'])
 'exchange-current density': 'single' (possible: ['single', 'current sigmoid'])
 'heat of mixing': 'false' (possible: ['false', 'true'])
 'hydrolysis': 'false' (possible: ['false', 'true'])
diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py
index 973a0f348b..8187729479 100644
--- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py
+++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_dfn.py
@@ -2,7 +2,6 @@
 # Tests for the lithium-ion DFN model
 #
 import pybamm
-import pytest
 from tests import BaseUnitTestLithiumIon
 
 
@@ -10,11 +9,6 @@ class TestDFN(BaseUnitTestLithiumIon):
     def setup_method(self):
         self.model = pybamm.lithium_ion.DFN
 
-    def test_electrolyte_options(self):
-        options = {"electrolyte conductivity": "integrated"}
-        with pytest.raises(pybamm.OptionError, match="electrolyte conductivity"):
-            pybamm.lithium_ion.DFN(options)
-
     def test_well_posed_size_distribution(self):
         options = {"particle size": "distribution"}
         self.check_well_posedness(options)
diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_newman_tobias.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_newman_tobias.py
index 58149a69de..6f35b35a82 100644
--- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_newman_tobias.py
+++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_newman_tobias.py
@@ -10,11 +10,6 @@ class TestNewmanTobias(BaseUnitTestLithiumIon):
     def setup_method(self):
         self.model = pybamm.lithium_ion.NewmanTobias
 
-    def test_electrolyte_options(self):
-        options = {"electrolyte conductivity": "integrated"}
-        with pytest.raises(pybamm.OptionError, match="electrolyte conductivity"):
-            pybamm.lithium_ion.NewmanTobias(options)
-
     @pytest.mark.skip(reason="Test currently not implemented")
     def test_well_posed_particle_phases(self):
         pass  # skip this test
diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_spme.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_spme.py
index ecab4384fc..38bde6c285 100644
--- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_spme.py
+++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/test_spme.py
@@ -32,7 +32,3 @@ def test_electrolyte_options(self):
         options = {"electrolyte conductivity": "full"}
         with pytest.raises(pybamm.OptionError, match="electrolyte conductivity"):
             pybamm.lithium_ion.SPMe(options)
-
-    def test_integrated_conductivity(self):
-        options = {"electrolyte conductivity": "integrated"}
-        self.check_well_posedness(options)