diff --git a/docs/source/models/submodels/active_material/index.rst b/docs/source/models/submodels/active_material/index.rst
index cff9324205..1022c9c6c0 100644
--- a/docs/source/models/submodels/active_material/index.rst
+++ b/docs/source/models/submodels/active_material/index.rst
@@ -8,6 +8,5 @@ Submodels for (loss of) active material
 
   base_active_material
   constant_active_material
-  reaction_driven_active_material
-  stress_driven_active_material
+  loss_active_material
 
diff --git a/docs/source/models/submodels/active_material/loss_active_material.rst b/docs/source/models/submodels/active_material/loss_active_material.rst
new file mode 100644
index 0000000000..b39a8d2448
--- /dev/null
+++ b/docs/source/models/submodels/active_material/loss_active_material.rst
@@ -0,0 +1,8 @@
+Loss of Active Material 
+=======================
+
+.. autoclass:: pybamm.active_material.LossActiveMaterial
+    :members:
+
+
+
diff --git a/docs/source/models/submodels/active_material/reaction_driven_active_material.rst b/docs/source/models/submodels/active_material/reaction_driven_active_material.rst
deleted file mode 100644
index d17f62ac6c..0000000000
--- a/docs/source/models/submodels/active_material/reaction_driven_active_material.rst
+++ /dev/null
@@ -1,8 +0,0 @@
-Reaction-driven Loss of Active Material
-=======================================
-
-.. autoclass:: pybamm.active_material.ReactionDriven
-    :members:
-
-
-
diff --git a/docs/source/models/submodels/active_material/stress_driven_active_material.rst b/docs/source/models/submodels/active_material/stress_driven_active_material.rst
deleted file mode 100644
index a4fa94fbd3..0000000000
--- a/docs/source/models/submodels/active_material/stress_driven_active_material.rst
+++ /dev/null
@@ -1,8 +0,0 @@
-Stress-driven Loss of Active Material 
-=====================================
-
-.. autoclass:: pybamm.active_material.StressDriven
-    :members:
-
-
-
diff --git a/pybamm/models/full_battery_models/base_battery_model.py b/pybamm/models/full_battery_models/base_battery_model.py
index c2f7a6f977..b88a732fca 100644
--- a/pybamm/models/full_battery_models/base_battery_model.py
+++ b/pybamm/models/full_battery_models/base_battery_model.py
@@ -58,7 +58,7 @@ class BatteryModelOptions(pybamm.FuzzyDict):
                 "reversible" or "irreversible".
             * "loss of active material" : str
                 Sets the model for loss of active material. Can be "none" (default),
-                "stress-driven", or "reaction-driven".
+                "stress-driven", "reaction-driven", or "stress and reaction-driven".
                 A 2-tuple can be provided for different behaviour in negative and
                 positive electrodes.
             * "operating mode" : str
@@ -174,7 +174,12 @@ def __init__(self, extra_options):
             "interface utilisation": ["full", "constant", "current-driven"],
             "lithium plating": ["none", "reversible", "irreversible"],
             "lithium plating porosity change": ["false", "true"],
-            "loss of active material": ["none", "stress-driven", "reaction-driven"],
+            "loss of active material": [
+                "none",
+                "stress-driven",
+                "reaction-driven",
+                "stress and reaction-driven",
+            ],
             "operating mode": ["current", "voltage", "power", "CCCV"],
             "particle": [
                 "Fickian diffusion",
@@ -235,7 +240,7 @@ def __init__(self, extra_options):
         # provided
         # return "none" if option not given
         lam_option = extra_options.get("loss of active material", "none")
-        if "stress-driven" in lam_option:
+        if "stress-driven" in lam_option or "stress and reaction-driven" in lam_option:
             default_options["particle mechanics"] = "swelling only"
         else:
             default_options["particle mechanics"] = "none"
@@ -423,6 +428,32 @@ def print_detailed_options(self):
         """
         print(self.__doc__)
 
+    @property
+    def negative(self):
+        "Returns the options for the negative electrode"
+        # index 0 in a 2-tuple for the negative electrode
+        return BatteryModelDomainOptions(self.items(), 0)
+
+    @property
+    def positive(self):
+        "Returns the options for the positive electrode"
+        # index 1 in a 2-tuple for the positive electrode
+        return BatteryModelDomainOptions(self.items(), 1)
+
+
+class BatteryModelDomainOptions(dict):
+    def __init__(self, dict_items, index):
+        super().__init__(dict_items)
+        self.index = index
+
+    def __getitem__(self, key):
+        options = super().__getitem__(key)
+        if isinstance(options, str):
+            return options
+        else:
+            # 2-tuple, first is negative domain, second is positive domain
+            return options[self.index]
+
 
 class BaseBatteryModel(pybamm.BaseModel):
     """
@@ -897,22 +928,15 @@ def set_current_collector_submodel(self):
         self.submodels["current collector"] = submodel
 
     def set_interface_utilisation_submodel(self):
-        # this option can either be a string (both sides the same) or a 2-tuple
-        # to indicate different options in negative and positive electrodes
-        if isinstance(self.options["interface utilisation"], str):
-            util_left = self.options["interface utilisation"]
-            util_right = self.options["interface utilisation"]
-        else:
-            util_left, util_right = self.options["interface utilisation"]
-
         if self.half_cell:
-            domains = [[util_left, "Counter"], [util_right, "Positive"]]
+            domains = ["Counter", "Positive"]
         else:
-            domains = [[util_left, "Negative"], [util_right, "Positive"]]
-        for util, domain in domains:
+            domains = ["Negative", "Positive"]
+        for domain in domains:
             name = domain.lower() + " interface utilisation"
             if domain == "Counter":
                 domain = "Negative"
+            util = getattr(self.options, domain.lower())["interface utilisation"]
             if util == "full":
                 self.submodels[name] = pybamm.interface_utilisation.Full(
                     self.param, domain, self.options
diff --git a/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py b/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py
index e6af73402a..0bf9b2661e 100644
--- a/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py
+++ b/pybamm/models/full_battery_models/lithium_ion/base_lithium_ion_model.py
@@ -225,21 +225,15 @@ def set_other_reaction_submodels_to_zero(self):
         )
 
     def set_crack_submodel(self):
-        # this option can either be a string (both sides the same) or a 2-tuple
-        # to indicate different options in negative and positive electrodes
-        if isinstance(self.options["particle mechanics"], str):
-            crack_left = self.options["particle mechanics"]
-            crack_right = self.options["particle mechanics"]
-        else:
-            crack_left, crack_right = self.options["particle mechanics"]
-        for crack_side, domain in [[crack_left, "Negative"], [crack_right, "Positive"]]:
-            if crack_side == "none":
+        for domain in ["Negative", "Positive"]:
+            crack = getattr(self.options, domain.lower())["particle mechanics"]
+            if crack == "none":
                 pass
-            elif crack_side == "swelling only":
+            elif crack == "swelling only":
                 self.submodels[
                     domain.lower() + " particle mechanics"
                 ] = pybamm.particle_mechanics.SwellingOnly(self.param, domain)
-            elif crack_side == "swelling and cracking":
+            elif crack == "swelling and cracking":
                 self.submodels[
                     domain.lower() + " particle mechanics"
                 ] = pybamm.particle_mechanics.CrackPropagation(
@@ -247,28 +241,16 @@ def set_crack_submodel(self):
                 )
 
     def set_active_material_submodel(self):
-        # this option can either be a string (both sides the same) or a 2-tuple
-        # to indicate different options in negative and positive electrodes
-        if isinstance(self.options["loss of active material"], str):
-            lam_left = self.options["loss of active material"]
-            lam_right = self.options["loss of active material"]
-        else:
-            lam_left, lam_right = self.options["loss of active material"]
-        for lam_side, domain in [[lam_left, "Negative"], [lam_right, "Positive"]]:
-            if lam_side == "none":
+        for domain in ["Negative", "Positive"]:
+            lam = getattr(self.options, domain.lower())["loss of active material"]
+            if lam == "none":
                 self.submodels[
                     domain.lower() + " active material"
                 ] = pybamm.active_material.Constant(self.param, domain, self.options)
-            elif lam_side == "stress-driven":
-                self.submodels[
-                    domain.lower() + " active material"
-                ] = pybamm.active_material.StressDriven(
-                    self.param, domain, self.options, self.x_average
-                )
-            elif lam_side == "reaction-driven":
+            else:
                 self.submodels[
                     domain.lower() + " active material"
-                ] = pybamm.active_material.ReactionDriven(
+                ] = pybamm.active_material.LossActiveMaterial(
                     self.param, domain, self.options, self.x_average
                 )
 
diff --git a/pybamm/models/full_battery_models/lithium_ion/dfn.py b/pybamm/models/full_battery_models/lithium_ion/dfn.py
index d92062dbff..8aa1866a6d 100644
--- a/pybamm/models/full_battery_models/lithium_ion/dfn.py
+++ b/pybamm/models/full_battery_models/lithium_ion/dfn.py
@@ -80,18 +80,10 @@ def set_intercalation_kinetics_submodel(self):
         )
 
     def set_particle_submodel(self):
-
-        if isinstance(self.options["particle"], str):
-            particle_left = self.options["particle"]
-            particle_right = self.options["particle"]
-        else:
-            particle_left, particle_right = self.options["particle"]
-        for particle_side, domain in [
-            [particle_left, "Negative"],
-            [particle_right, "Positive"],
-        ]:
+        for domain in ["Negative", "Positive"]:
+            particle = getattr(self.options, domain.lower())["particle"]
             if self.options["particle size"] == "single":
-                if particle_side == "Fickian diffusion":
+                if particle == "Fickian diffusion":
                     self.submodels[
                         domain.lower() + " particle"
                     ] = pybamm.particle.no_distribution.FickianDiffusion(
@@ -99,7 +91,7 @@ def set_particle_submodel(self):
                         domain,
                         self.options,
                     )
-                elif particle_side in [
+                elif particle in [
                     "uniform profile",
                     "quadratic profile",
                     "quartic profile",
@@ -107,19 +99,16 @@ def set_particle_submodel(self):
                     self.submodels[
                         domain.lower() + " particle"
                     ] = pybamm.particle.no_distribution.PolynomialProfile(
-                        self.param,
-                        domain,
-                        particle_side,
-                        self.options,
+                        self.param, domain, particle, self.options
                     )
             elif self.options["particle size"] == "distribution":
-                if particle_side == "Fickian diffusion":
+                if particle == "Fickian diffusion":
                     self.submodels[
                         domain.lower() + " particle"
                     ] = pybamm.particle.size_distribution.FickianDiffusion(
                         self.param, domain
                     )
-                elif particle_side == "uniform profile":
+                elif particle == "uniform profile":
                     self.submodels[
                         domain.lower() + " particle"
                     ] = pybamm.particle.size_distribution.UniformProfile(
diff --git a/pybamm/models/full_battery_models/lithium_ion/newman_tobias.py b/pybamm/models/full_battery_models/lithium_ion/newman_tobias.py
index 877f8de420..ba3be1c9be 100644
--- a/pybamm/models/full_battery_models/lithium_ion/newman_tobias.py
+++ b/pybamm/models/full_battery_models/lithium_ion/newman_tobias.py
@@ -56,22 +56,15 @@ def __init__(self, options=None, name="Newman-Tobias model", build=True):
         pybamm.citations.register("Chu2020")
 
     def set_particle_submodel(self):
-        if isinstance(self.options["particle"], str):
-            particle_left = self.options["particle"]
-            particle_right = self.options["particle"]
-        else:
-            particle_left, particle_right = self.options["particle"]
-        for particle_side, domain in [
-            [particle_left, "Negative"],
-            [particle_right, "Positive"],
-        ]:
-            if particle_side == "Fickian diffusion":
+        for domain in ["Negative", "Positive"]:
+            particle = getattr(self.options, domain.lower())["particle"]
+            if particle == "Fickian diffusion":
                 self.submodels[
                     domain.lower() + " particle"
                 ] = pybamm.particle.no_distribution.XAveragedFickianDiffusion(
                     self.param, domain, self.options
                 )
-            elif particle_side in [
+            elif particle in [
                 "uniform profile",
                 "quadratic profile",
                 "quartic profile",
@@ -79,7 +72,7 @@ def set_particle_submodel(self):
                 self.submodels[
                     domain.lower() + " particle"
                 ] = pybamm.particle.no_distribution.XAveragedPolynomialProfile(
-                    self.param, domain, particle_side, self.options
+                    self.param, domain, particle, self.options
                 )
 
     def set_electrolyte_submodel(self):
diff --git a/pybamm/models/full_battery_models/lithium_ion/spm.py b/pybamm/models/full_battery_models/lithium_ion/spm.py
index 44a95dc690..73d1a26347 100644
--- a/pybamm/models/full_battery_models/lithium_ion/spm.py
+++ b/pybamm/models/full_battery_models/lithium_ion/spm.py
@@ -107,22 +107,15 @@ def set_intercalation_kinetics_submodel(self):
             )
 
     def set_particle_submodel(self):
-        if isinstance(self.options["particle"], str):
-            particle_left = self.options["particle"]
-            particle_right = self.options["particle"]
-        else:
-            particle_left, particle_right = self.options["particle"]
-        for particle_side, domain in [
-            [particle_left, "Negative"],
-            [particle_right, "Positive"],
-        ]:
-            if particle_side == "Fickian diffusion":
+        for domain in ["Negative", "Positive"]:
+            particle = getattr(self.options, domain.lower())["particle"]
+            if particle == "Fickian diffusion":
                 self.submodels[
                     domain.lower() + " particle"
                 ] = pybamm.particle.no_distribution.XAveragedFickianDiffusion(
                     self.param, domain, self.options
                 )
-            elif particle_side in [
+            elif particle in [
                 "uniform profile",
                 "quadratic profile",
                 "quartic profile",
@@ -130,7 +123,7 @@ def set_particle_submodel(self):
                 self.submodels[
                     domain.lower() + " particle"
                 ] = pybamm.particle.no_distribution.XAveragedPolynomialProfile(
-                    self.param, domain, particle_side, self.options
+                    self.param, domain, particle, self.options
                 )
 
     def set_solid_submodel(self):
diff --git a/pybamm/models/submodels/active_material/__init__.py b/pybamm/models/submodels/active_material/__init__.py
index c8031a1a6a..1fd47a39d5 100644
--- a/pybamm/models/submodels/active_material/__init__.py
+++ b/pybamm/models/submodels/active_material/__init__.py
@@ -1,4 +1,3 @@
 from .base_active_material import BaseModel
 from .constant_active_material import Constant
-from .stress_driven_active_material import StressDriven
-from .reaction_driven_active_material import ReactionDriven
+from .loss_active_material import LossActiveMaterial
diff --git a/pybamm/models/submodels/active_material/stress_driven_active_material.py b/pybamm/models/submodels/active_material/loss_active_material.py
similarity index 55%
rename from pybamm/models/submodels/active_material/stress_driven_active_material.py
rename to pybamm/models/submodels/active_material/loss_active_material.py
index a477539061..4cdccbc414 100644
--- a/pybamm/models/submodels/active_material/stress_driven_active_material.py
+++ b/pybamm/models/submodels/active_material/loss_active_material.py
@@ -6,9 +6,8 @@
 from .base_active_material import BaseModel
 
 
-class StressDriven(BaseModel):
-    """Submodel for varying active material volume fraction, driven by stress, from
-    [1]_ and [2]_.
+class LossActiveMaterial(BaseModel):
+    """Submodel for varying active material volume fraction from [1]_ and [2]_.
 
     Parameters
     ----------
@@ -56,42 +55,76 @@ def get_fundamental_variables(self):
         return variables
 
     def get_coupled_variables(self, variables):
-        # obtain the rate of loss of active materials (LAM) by stress
-        # This is loss of active material model by mechanical effects
-        if self.x_average is True:
-            stress_t_surf = variables[
-                "X-averaged "
-                + self.domain.lower()
-                + " particle surface tangential stress"
-            ]
-            stress_r_surf = variables[
-                "X-averaged " + self.domain.lower() + " particle surface radial stress"
-            ]
-        else:
-            stress_t_surf = variables[
-                self.domain + " particle surface tangential stress"
-            ]
-            stress_r_surf = variables[self.domain + " particle surface radial stress"]
-
-        if self.domain == "Negative":
-            beta_LAM = self.param.beta_LAM_n
-            stress_critical = self.param.stress_critical_n
-            m_LAM = self.param.m_LAM_n
-        else:
-            beta_LAM = self.param.beta_LAM_p
-            stress_critical = self.param.stress_critical_p
-            m_LAM = self.param.m_LAM_p
-
-        stress_h_surf = (stress_r_surf + 2 * stress_t_surf) / 3
-        # compressive stress make no contribution
-        stress_h_surf *= stress_h_surf > 0
-        # assuming the minimum hydrostatic stress is zero for full cycles
-        stress_h_surf_min = stress_h_surf * 0
-        j_stress_LAM = (
-            -beta_LAM * ((stress_h_surf - stress_h_surf_min) / stress_critical) ** m_LAM
-        )
-
-        deps_solid_dt = j_stress_LAM
+        deps_solid_dt = 0
+        lam_option = getattr(self.options, self.domain.lower())[
+            "loss of active material"
+        ]
+        if "stress" in lam_option:
+            # obtain the rate of loss of active materials (LAM) by stress
+            # This is loss of active material model by mechanical effects
+            if self.x_average is True:
+                stress_t_surf = variables[
+                    "X-averaged "
+                    + self.domain.lower()
+                    + " particle surface tangential stress"
+                ]
+                stress_r_surf = variables[
+                    "X-averaged "
+                    + self.domain.lower()
+                    + " particle surface radial stress"
+                ]
+            else:
+                stress_t_surf = variables[
+                    self.domain + " particle surface tangential stress"
+                ]
+                stress_r_surf = variables[
+                    self.domain + " particle surface radial stress"
+                ]
+
+            if self.domain == "Negative":
+                beta_LAM = self.param.beta_LAM_n
+                stress_critical = self.param.stress_critical_n
+                m_LAM = self.param.m_LAM_n
+            else:
+                beta_LAM = self.param.beta_LAM_p
+                stress_critical = self.param.stress_critical_p
+                m_LAM = self.param.m_LAM_p
+
+            stress_h_surf = (stress_r_surf + 2 * stress_t_surf) / 3
+            # compressive stress make no contribution
+            stress_h_surf *= stress_h_surf > 0
+            # assuming the minimum hydrostatic stress is zero for full cycles
+            stress_h_surf_min = stress_h_surf * 0
+            j_stress_LAM = (
+                -beta_LAM
+                * ((stress_h_surf - stress_h_surf_min) / stress_critical) ** m_LAM
+            )
+            deps_solid_dt += j_stress_LAM
+
+        if "reaction" in lam_option:
+            if self.x_average is True:
+                a = variables[
+                    "X-averaged "
+                    + self.domain.lower()
+                    + " electrode surface area to volume ratio"
+                ]
+            else:
+                a = variables[self.domain + " electrode surface area to volume ratio"]
+
+            if self.domain == "Negative":
+                beta_LAM_sei = self.param.beta_LAM_sei_n
+                if self.x_average is True:
+                    j_sei = variables["X-averaged SEI interfacial current density"]
+                else:
+                    j_sei = variables["SEI interfacial current density"]
+            else:
+                # No SEI in the positive electrode so no reaction-driven LAM
+                # until other reactions are implemented
+                beta_LAM_sei = self.param.beta_LAM_sei_p
+                j_sei = 0
+
+            j_stress_reaction = beta_LAM_sei * a * j_sei
+            deps_solid_dt += j_stress_reaction
         variables.update(
             self._get_standard_active_material_change_variables(deps_solid_dt)
         )
diff --git a/pybamm/models/submodels/active_material/reaction_driven_active_material.py b/pybamm/models/submodels/active_material/reaction_driven_active_material.py
deleted file mode 100644
index 0a1a1be844..0000000000
--- a/pybamm/models/submodels/active_material/reaction_driven_active_material.py
+++ /dev/null
@@ -1,120 +0,0 @@
-#
-# Class for varying active material volume fraction, driven by reactions
-#
-import pybamm
-
-from .base_active_material import BaseModel
-
-
-class ReactionDriven(BaseModel):
-    """Submodel for varying active material volume fraction, driven by reactions, from
-    [1]_
-
-    Parameters
-    ----------
-    param : parameter class
-        The parameters to use for this submodel
-    domain : str
-        The domain of the model either 'Negative' or 'Positive'
-    options : dict
-        Additional options to pass to the model
-    x_average : bool
-        Whether to use x-averaged variables (SPM, SPMe, etc) or full variables (DFN)
-
-    **Extends:** :class:`pybamm.active_material.BaseModel`
-
-    References
-    ----------
-    .. [1] Reniers, J. M., Mulder, G., & Howey, D. A. (2019). Review and performance
-           comparison of mechanical-chemical degradation models for lithium-ion
-           batteries. Journal of The Electrochemical Society, 166(14), A3189.
-    """
-
-    def __init__(self, param, domain, options, x_average):
-        super().__init__(param, domain, options=options)
-        pybamm.citations.register("Reniers2019")
-        self.x_average = x_average
-
-    def get_fundamental_variables(self):
-        domain = self.domain.lower() + " electrode"
-        if self.x_average is True:
-            eps_solid_xav = pybamm.Variable(
-                "X-averaged " + domain + " active material volume fraction",
-                domain="current collector",
-            )
-            eps_solid = pybamm.PrimaryBroadcast(eps_solid_xav, domain)
-        else:
-            eps_solid = pybamm.Variable(
-                self.domain + " electrode active material volume fraction",
-                domain=domain,
-                auxiliary_domains={"secondary": "current collector"},
-            )
-        variables = self._get_standard_active_material_variables(eps_solid)
-        return variables
-
-    def get_coupled_variables(self, variables):
-        if self.x_average is True:
-            a = variables[
-                "X-averaged "
-                + self.domain.lower()
-                + " electrode surface area to volume ratio"
-            ]
-        else:
-            a = variables[self.domain + " electrode surface area to volume ratio"]
-
-        if self.domain == "Negative":
-            beta_LAM_sei = self.param.beta_LAM_sei_n
-            if self.x_average is True:
-                j_sei = variables["X-averaged SEI interfacial current density"]
-            else:
-                j_sei = variables["SEI interfacial current density"]
-        else:
-            beta_LAM_sei = self.param.beta_LAM_sei_p
-            j_sei = 0
-
-        deps_solid_dt = beta_LAM_sei * a * j_sei
-        variables.update(
-            self._get_standard_active_material_change_variables(deps_solid_dt)
-        )
-        return variables
-
-    def set_rhs(self, variables):
-        domain = self.domain.lower() + " electrode"
-        if self.x_average is True:
-            eps_solid = variables[
-                "X-averaged " + domain + " active material volume fraction"
-            ]
-            deps_solid_dt = variables[
-                "X-averaged " + domain + " active material volume fraction change"
-            ]
-        else:
-            eps_solid = variables[
-                self.domain + " electrode active material volume fraction"
-            ]
-            deps_solid_dt = variables[
-                self.domain + " electrode active material volume fraction change"
-            ]
-
-        self.rhs = {eps_solid: deps_solid_dt}
-
-    def set_initial_conditions(self, variables):
-
-        if self.domain == "Negative":
-            x_n = pybamm.standard_spatial_vars.x_n
-            eps_solid_init = self.param.epsilon_s_n(x_n)
-        elif self.domain == "Positive":
-            x_p = pybamm.standard_spatial_vars.x_p
-            eps_solid_init = self.param.epsilon_s_p(x_p)
-
-        if self.x_average is True:
-            eps_solid_xav = variables[
-                "X-averaged "
-                + self.domain.lower()
-                + " electrode active material volume fraction"
-            ]
-            self.initial_conditions = {eps_solid_xav: pybamm.x_average(eps_solid_init)}
-        else:
-            eps_solid = variables[
-                self.domain + " electrode active material volume fraction"
-            ]
-            self.initial_conditions = {eps_solid: eps_solid_init}
diff --git a/pybamm/models/submodels/particle/no_distribution/base_fickian.py b/pybamm/models/submodels/particle/no_distribution/base_fickian.py
index b61be85394..3e2634d421 100644
--- a/pybamm/models/submodels/particle/no_distribution/base_fickian.py
+++ b/pybamm/models/submodels/particle/no_distribution/base_fickian.py
@@ -36,22 +36,19 @@ def _get_effective_diffusivity(self, c, T):
 
         # Account for stress-induced diffusion by defining a multiplicative
         # "stress factor"
-        # This option can either be a string (both sides the same) or a 2-tuple
-        # to indicate different options in negative and positive electrodes
-        if isinstance(self.options["stress-induced diffusion"], str):
-            stress_left = self.options["stress-induced diffusion"]
-            stress_right = self.options["stress-induced diffusion"]
-        else:
-            stress_left, stress_right = self.options["stress-induced diffusion"]
+        stress_option = getattr(self.options, self.domain.lower())[
+            "stress-induced diffusion"
+        ]
 
-        if self.domain == "Negative" and stress_left == "true":
-            stress_factor = 1 + param.theta_n * (c - param.c_n_0) / (
-                1 + param.Theta * T
-            )
-        elif self.domain == "Positive" and stress_right == "true":
-            stress_factor = 1 + param.theta_p * (c - param.c_p_0) / (
-                1 + param.Theta * T
-            )
+        if stress_option == "true":
+            if self.domain == "Negative":
+                stress_factor = 1 + param.theta_n * (c - param.c_n_0) / (
+                    1 + param.Theta * T
+                )
+            elif self.domain == "Positive":
+                stress_factor = 1 + param.theta_p * (c - param.c_p_0) / (
+                    1 + param.Theta * T
+                )
         else:
             stress_factor = 1
 
diff --git a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py
index 0a2558fc3f..ac09973c5c 100644
--- a/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py
+++ b/tests/integration/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py
@@ -121,10 +121,6 @@ def test_current_driven_utilisation(self):
         )
         self.run_basic_processing_test(options, parameter_values=parameter_values)
 
-    def test_loss_active_material_reaction_both(self):
-        options = {"loss of active material": "reaction-driven"}
-        self.run_basic_processing_test(options)
-
     def test_surface_form_differential(self):
         options = {"surface form": "differential"}
         self.run_basic_processing_test(options)
@@ -205,6 +201,15 @@ def test_loss_active_material_stress_both(self):
         parameter_values = pybamm.ParameterValues("Ai2020")
         self.run_basic_processing_test(options, parameter_values=parameter_values)
 
+    def test_loss_active_material_reaction(self):
+        options = {"loss of active material": "reaction-driven"}
+        self.run_basic_processing_test(options)
+
+    def test_loss_active_material_stress_and_reaction(self):
+        options = {"loss of active material": "stress and reaction-driven"}
+        parameter_values = pybamm.ParameterValues("Ai2020")
+        self.run_basic_processing_test(options, parameter_values=parameter_values)
+
     def test_negative_cracking(self):
         options = {"particle mechanics": ("swelling and cracking", "none")}
         parameter_values = pybamm.ParameterValues("Ai2020")
diff --git a/tests/unit/test_citations.py b/tests/unit/test_citations.py
index 6e9358895e..0861479099 100644
--- a/tests/unit/test_citations.py
+++ b/tests/unit/test_citations.py
@@ -168,12 +168,7 @@ def test_reniers_2019(self):
 
         citations._reset()
         self.assertNotIn("Reniers2019", citations._papers_to_cite)
-        pybamm.active_material.StressDriven(None, None, None, True)
-        self.assertIn("Reniers2019", citations._papers_to_cite)
-
-        citations._reset()
-        self.assertNotIn("Reniers2019", citations._papers_to_cite)
-        pybamm.active_material.ReactionDriven(None, None, None, True)
+        pybamm.active_material.LossActiveMaterial(None, None, None, True)
         self.assertIn("Reniers2019", citations._papers_to_cite)
 
     def test_mohtat_2019(self):
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 6882ea458e..cd452ebdd2 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 @@
 'interface utilisation': 'full' (possible: ['full', 'constant', 'current-driven'])
 'lithium plating': 'none' (possible: ['none', 'reversible', 'irreversible'])
 'lithium plating porosity change': 'false' (possible: ['false', 'true'])
-'loss of active material': 'stress-driven' (possible: ['none', 'stress-driven', 'reaction-driven'])
+'loss of active material': 'stress-driven' (possible: ['none', 'stress-driven', 'reaction-driven', 'stress and reaction-driven'])
 'operating mode': 'current' (possible: ['current', 'voltage', 'power', 'CCCV'])
 'particle': 'Fickian diffusion' (possible: ['Fickian diffusion', 'fast diffusion', 'uniform profile', 'quadratic profile', 'quartic profile'])
 'particle mechanics': 'swelling only' (possible: ['none', 'swelling only', 'swelling and cracking'])
@@ -322,6 +322,16 @@ def test_print_options(self):
             output = buffer.getvalue()
         self.assertEqual(output, PRINT_OPTIONS_OUTPUT)
 
+    def test_domain_options(self):
+        options = BatteryModelOptions(
+            {"particle": ("Fickian diffusion", "quadratic profile")}
+        )
+        self.assertEqual(options.negative["particle"], "Fickian diffusion")
+        self.assertEqual(options.positive["particle"], "quadratic profile")
+        # something that is the same in both domains
+        self.assertEqual(options.negative["thermal"], "isothermal")
+        self.assertEqual(options.positive["thermal"], "isothermal")
+
 
 if __name__ == "__main__":
     print("Add -v for more debug output")
diff --git a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py
index ddbcf55806..346ae9e8ff 100644
--- a/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py
+++ b/tests/unit/test_models/test_full_battery_models/test_lithium_ion/base_lithium_ion_tests.py
@@ -100,10 +100,14 @@ def test_well_posed_loss_active_material_stress_both(self):
         options = {"loss of active material": "stress-driven"}
         self.check_well_posedness(options)
 
-    def test_well_posed_loss_active_material_stress_reaction_both(self):
+    def test_well_posed_loss_active_material_reaction(self):
         options = {"loss of active material": "reaction-driven"}
         self.check_well_posedness(options)
 
+    def test_well_posed_loss_active_material_stress_reaction(self):
+        options = {"loss of active material": "stress and reaction-driven"}
+        self.check_well_posedness(options)
+
     def test_well_posed_surface_form_differential(self):
         options = {"surface form": "differential"}
         self.check_well_posedness(options)