Repository: pybamm-team/PyBaMM

Issue:
Improve memory usage for experiments with `start_time`
At the moment, two steps with the same operating conditions but different `start_time` are treated as different. This means that the model needs to be processed twice, which can end up with a big memory consumption (and longer running times) for very long experiments.

A way to improve this would be to read the `start_time` from the unprocessed step, this way we can have a single processed step for all of them.

Discussion:
review-notebook-app[bot]: Check out this pull request on&nbsp; <a href="https://app.reviewnb.com/pybamm-team/PyBaMM/pull/3261"><img align="absmiddle"  alt="ReviewNB" height="28" class="BotMessageButtonImage" src="https://raw.githubusercontent.com/ReviewNB/support/master/images/button_reviewnb.png"/></a> 

 See visual diffs & provide feedback on Jupyter Notebooks. 

---

 <i>Powered by <a href='https://www.reviewnb.com/?utm_source=gh'>ReviewNB</a></i>

Patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index b4f13b48a..971b248ab 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -11,6 +11,7 @@
 - Fixed a bug that caused incorrect results of “{Domain} electrode thickness change [m]” due to the absence of dimension for the variable `electrode_thickness_change`([#3329](https://github.com/pybamm-team/PyBaMM/pull/3329)).
 - Fixed a bug that occured in `check_ys_are_not_too_large` when trying to reference `y-slice` where the referenced variable was not a `pybamm.StateVector` ([#3313](https://github.com/pybamm-team/PyBaMM/pull/3313)
 - Fixed a bug with `_Heaviside._evaluate_for_shape` which meant some expressions involving heaviside function and subtractions did not work ([#3306](https://github.com/pybamm-team/PyBaMM/pull/3306))
+- Fixed bug causing incorrect activation energies using `create_from_bpx()` ([#3242](https://github.com/pybamm-team/PyBaMM/pull/3242))
 - The `OneDimensionalX` thermal model has been updated to account for edge/tab cooling and account for the current collector volumetric heat capacity. It now gives the correct behaviour compared with a lumped model with the correct total heat transfer coefficient and surface area for cooling. ([#3042](https://github.com/pybamm-team/PyBaMM/pull/3042))
 - Fixed a bug where the "basic" lithium-ion models gave incorrect results when using nonlinear particle diffusivity ([#3207](https://github.com/pybamm-team/PyBaMM/pull/3207))
 - Particle size distributions now work with SPMe and NewmanTobias models ([#3207](https://github.com/pybamm-team/PyBaMM/pull/3207))
@@ -22,12 +23,15 @@
 - Error generated when invalid parameter values are passed ([#3132](https://github.com/pybamm-team/PyBaMM/pull/3132))
 - Parameters in `Prada2013` have been updated to better match those given in the paper, which is a 2.3 Ah cell, instead of the mix-and-match with the 1.1 Ah cell from Lain2019 ([#3096](https://github.com/pybamm-team/PyBaMM/pull/3096))
 
+## Optimizations
 
+- Improved how steps are processed in simulations to reduce memory usage ([#3261](https://github.com/pybamm-team/PyBaMM/pull/3261))
 
 ## Breaking changes
 
 - The class `pybamm.thermal.OneDimensionalX` has been moved to `pybamm.thermal.pouch_cell.OneDimensionalX` to reflect the fact that the model formulation implicitly assumes a pouch cell geometry ([#3257](https://github.com/pybamm-team/PyBaMM/pull/3257))
 - The "lumped" thermal option now always used the parameters "Cell cooling surface area [m2]", "Cell volume [m3]" and "Total heat transfer coefficient [W.m-2.K-1]" to compute the cell cooling regardless of the chosen "cell geometry" option. The user must now specify the correct values for these parameters instead of them being calculated based on e.g. a pouch cell. An `OptionWarning` is raised to let users know to update their parameters ([#3257](https://github.com/pybamm-team/PyBaMM/pull/3257))
+- Numpy functions now work with PyBaMM symbols (e.g. `np.exp(pybamm.Symbol("a"))` returns `pybamm.Exp(pybamm.Symbol("a"))`). This means that parameter functions can be specified using numpy functions instead of pybamm functions. Additionally, combining numpy arrays with pybamm objects now works (the numpy array is converted to a pybamm array) ([#3205](https://github.com/pybamm-team/PyBaMM/pull/3205))
 - Added option to use an empirical hysteresis model for the diffusivity and exchange-current density ([#3194](https://github.com/pybamm-team/PyBaMM/pull/3194))
 - Double-layer capacity can now be provided as a function of temperature ([#3174](https://github.com/pybamm-team/PyBaMM/pull/3174))
 - `pybamm_install_jax` is deprecated. It is now replaced with `pip install pybamm[jax]` ([#3163](https://github.com/pybamm-team/PyBaMM/pull/3163))
diff --git a/pybamm/experiment/experiment.py b/pybamm/experiment/experiment.py
index edcaeb8f5..d1c45015b 100644
--- a/pybamm/experiment/experiment.py
+++ b/pybamm/experiment/experiment.py
@@ -68,11 +68,6 @@ class Experiment:
             termination,
         )
 
-        self.datetime_formats = [
-            "Day %j %H:%M:%S",
-            "%Y-%m-%d %H:%M:%S",
-        ]
-
         operating_conditions_cycles = []
         for cycle in operating_conditions:
             # Check types and convert to list
@@ -89,21 +84,36 @@ class Experiment:
 
         # Convert strings to pybamm.step._Step objects
         # We only do this once per unique step, do avoid unnecessary conversions
-        unique_steps_unprocessed = set(operating_conditions_steps_unprocessed)
+        # Assign experiment period and temperature if not specified in step
+        self.period = _convert_time_to_seconds(period)
+        self.temperature = _convert_temperature_to_kelvin(temperature)
+
         processed_steps = {}
-        for step in unique_steps_unprocessed:
-            if isinstance(step, str):
-                processed_steps[step] = pybamm.step.string(step)
+        for step in operating_conditions_steps_unprocessed:
+            if repr(step) in processed_steps:
+                continue
+            elif isinstance(step, str):
+                processed_step = pybamm.step.string(step)
             elif isinstance(step, pybamm.step._Step):
-                processed_steps[step] = step
+                processed_step = step
+
+            if processed_step.period is None:
+                processed_step.period = self.period
+            if processed_step.temperature is None:
+                processed_step.temperature = self.temperature
+
+            processed_steps[repr(step)] = processed_step
+
+        self.operating_conditions_steps = [
+            processed_steps[repr(step)]
+            for step in operating_conditions_steps_unprocessed
+        ]
 
         # Save the processed unique steps and the processed operating conditions
         # for every step
         self.unique_steps = set(processed_steps.values())
-        self.operating_conditions_steps = [
-            processed_steps[step] for step in operating_conditions_steps_unprocessed
-        ]
 
+        # Allocate experiment global variables
         self.initial_start_time = self.operating_conditions_steps[0].start_time
 
         if (
@@ -118,15 +128,6 @@ class Experiment:
         self.termination_string = termination
         self.termination = self.read_termination(termination)
 
-        # Modify steps with period and temperature in place
-        self.period = _convert_time_to_seconds(period)
-        self.temperature = _convert_temperature_to_kelvin(temperature)
-        for step in self.unique_steps:
-            if step.period is None:
-                step.period = self.period
-            if step.temperature is None:
-                step.temperature = self.temperature
-
     def __str__(self):
         return str(self.operating_conditions_cycles)
 
diff --git a/pybamm/parameters/bpx.py b/pybamm/parameters/bpx.py
index 8f555fa9a..8efd26cd5 100644
--- a/pybamm/parameters/bpx.py
+++ b/pybamm/parameters/bpx.py
@@ -119,6 +119,9 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
     # reference temperature
     T_ref = pybamm_dict["Reference temperature [K]"]
 
+    def arrhenius(Ea, T):
+        return exp(Ea / constants.R * (1 / T_ref - 1 / T))
+
     # lumped parameters
     for name in [
         "Specific heat capacity [J.K-1.kg-1]",
@@ -258,10 +261,12 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
         "Maximum concentration in " + negative_electrode.pre_name.lower() + "[mol.m-3]"
     ]
     k_n_norm = pybamm_dict[
-        negative_electrode.pre_name + "reaction rate constant [mol.m-2.s-1]"
+        negative_electrode.pre_name
+        + "reaction rate constant [mol.m-2.s-1]"
     ]
-    E_a_n = pybamm_dict.get(
-        negative_electrode.pre_name + "reaction rate activation energy [J.mol-1]", 0.0
+    Ea_k_n = pybamm_dict.get(
+        negative_electrode.pre_name
+        + "reaction rate constant activation energy [J.mol-1]", 0.0
     )
     # Note that in BPX j = 2*F*k_norm*sqrt((ce/ce0)*(c/c_max)*(1-c/c_max))*sinh(...),
     # and in PyBaMM j = 2*k*sqrt(ce*c*(c_max - c))*sinh(...)
@@ -270,10 +275,9 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
     def _negative_electrode_exchange_current_density(c_e, c_s_surf, c_s_max, T):
         k_ref = k_n  # (A/m2)(m3/mol)**1.5 - includes ref concentrations
 
-        arrhenius = exp(E_a_n / constants.R * (1 / T_ref - 1 / T))
         return (
             k_ref
-            * arrhenius
+            * arrhenius(Ea_k_n, T)
             * c_e**0.5
             * c_s_surf**0.5
             * (c_s_max - c_s_surf) ** 0.5
@@ -288,10 +292,12 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
         "Maximum concentration in " + positive_electrode.pre_name.lower() + "[mol.m-3]"
     ]
     k_p_norm = pybamm_dict[
-        positive_electrode.pre_name + "reaction rate constant [mol.m-2.s-1]"
+        positive_electrode.pre_name
+        + "reaction rate constant [mol.m-2.s-1]"
     ]
-    E_a_p = pybamm_dict.get(
-        positive_electrode.pre_name + "reaction rate activation energy [J.mol-1]", 0.0
+    Ea_k_p = pybamm_dict.get(
+        positive_electrode.pre_name
+        + "reaction rate constant activation energy [J.mol-1]", 0.0
     )
     # Note that in BPX j = 2*F*k_norm*sqrt((ce/ce0)*(c/c_max)*(1-c/c_max))*sinh(...),
     # and in PyBaMM j = 2*k*sqrt(ce*c*(c_max - c))*sinh(...)
@@ -300,10 +306,9 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
     def _positive_electrode_exchange_current_density(c_e, c_s_surf, c_s_max, T):
         k_ref = k_p  # (A/m2)(m3/mol)**1.5 - includes ref concentrations
 
-        arrhenius = exp(E_a_p / constants.R * (1 / T_ref - 1 / T))
         return (
             k_ref
-            * arrhenius
+            * arrhenius(Ea_k_p, T)
             * c_e**0.5
             * c_s_surf**0.5
             * (c_s_max - c_s_surf) ** 0.5
@@ -316,7 +321,7 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
     # diffusivity
 
     # negative electrode
-    E_a = pybamm_dict.get(
+    Ea_D_n = pybamm_dict.get(
         negative_electrode.pre_name + "diffusivity activation energy [J.mol-1]", 0.0
     )
     D_n_ref = pybamm_dict[negative_electrode.pre_name + "diffusivity [m2.s-1]"]
@@ -324,30 +329,27 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
     if callable(D_n_ref):
 
         def _negative_electrode_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * D_n_ref(sto)
+            return arrhenius(Ea_D_n, T) * D_n_ref(sto)
 
     elif isinstance(D_n_ref, tuple):
 
         def _negative_electrode_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
             name, (x, y) = D_n_ref
-            return arrhenius * pybamm.Interpolant(
+            return arrhenius(Ea_D_n, T) * pybamm.Interpolant(
                 x, y, sto, name=name, interpolator="linear"
             )
 
     else:
 
         def _negative_electrode_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * D_n_ref
+            return arrhenius(Ea_D_n, T) * D_n_ref
 
     pybamm_dict[negative_electrode.pre_name + "diffusivity [m2.s-1]"] = _copy_func(
         _negative_electrode_diffusivity
     )
 
     # positive electrode
-    E_a = pybamm_dict.get(
+    Ea_D_p = pybamm_dict.get(
         positive_electrode.pre_name + "diffusivity activation energy [J.mol-1]", 0.0
     )
     D_p_ref = pybamm_dict[positive_electrode.pre_name + "diffusivity [m2.s-1]"]
@@ -355,30 +357,27 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
     if callable(D_p_ref):
 
         def _positive_electrode_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * D_p_ref(sto)
+            return arrhenius(Ea_D_p, T) * D_p_ref(sto)
 
     elif isinstance(D_p_ref, tuple):
 
         def _positive_electrode_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
             name, (x, y) = D_p_ref
-            return arrhenius * pybamm.Interpolant(
+            return arrhenius(Ea_D_p, T) * pybamm.Interpolant(
                 x, y, sto, name=name, interpolator="linear"
             )
 
     else:
 
         def _positive_electrode_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * D_p_ref
+            return arrhenius(Ea_D_p, T) * D_p_ref
 
     pybamm_dict[positive_electrode.pre_name + "diffusivity [m2.s-1]"] = _copy_func(
         _positive_electrode_diffusivity
     )
 
     # electrolyte
-    E_a = pybamm_dict.get(
+    Ea_D_e = pybamm_dict.get(
         electrolyte.pre_name + "diffusivity activation energy [J.mol-1]", 0.0
     )
     D_e_ref = pybamm_dict[electrolyte.pre_name + "diffusivity [m2.s-1]"]
@@ -386,54 +385,48 @@ def _bpx_to_param_dict(bpx: BPX) -> dict:
     if callable(D_e_ref):
 
         def _electrolyte_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * D_e_ref(sto)
+            return arrhenius(Ea_D_e, T) * D_e_ref(sto)
 
     elif isinstance(D_e_ref, tuple):
 
         def _electrolyte_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
             name, (x, y) = D_e_ref
-            return arrhenius * pybamm.Interpolant(
+            return arrhenius(Ea_D_e, T) * pybamm.Interpolant(
                 x, y, sto, name=name, interpolator="linear"
             )
 
     else:
 
         def _electrolyte_diffusivity(sto, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * D_e_ref
+            return arrhenius(Ea_D_e, T) * D_e_ref
 
     pybamm_dict[electrolyte.pre_name + "diffusivity [m2.s-1]"] = _copy_func(
         _electrolyte_diffusivity
     )
 
     # conductivity
-    E_a = pybamm_dict.get(
+    Ea_sigma_e = pybamm_dict.get(
         electrolyte.pre_name + "conductivity activation energy [J.mol-1]", 0.0
     )
-    C_e_ref = pybamm_dict[electrolyte.pre_name + "conductivity [S.m-1]"]
+    sigma_e_ref = pybamm_dict[electrolyte.pre_name + "conductivity [S.m-1]"]
 
-    if callable(C_e_ref):
+    if callable(sigma_e_ref):
 
         def _conductivity(c_e, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * C_e_ref(c_e)
+            return arrhenius(Ea_sigma_e, T) * sigma_e_ref(c_e)
 
-    elif isinstance(C_e_ref, tuple):
+    elif isinstance(sigma_e_ref, tuple):
 
         def _conductivity(c_e, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            name, (x, y) = C_e_ref
-            return arrhenius * pybamm.Interpolant(
+            name, (x, y) = sigma_e_ref
+            return arrhenius(Ea_sigma_e, T) * pybamm.Interpolant(
                 x, y, c_e, name=name, interpolator="linear"
             )
 
     else:
 
         def _conductivity(c_e, T):
-            arrhenius = exp(E_a / constants.R * (1 / T_ref - 1 / T))
-            return arrhenius * C_e_ref
+            return arrhenius(Ea_sigma_e, T) * sigma_e_ref
 
     pybamm_dict[electrolyte.pre_name + "conductivity [S.m-1]"] = _copy_func(
         _conductivity
diff --git a/pybamm/simulation.py b/pybamm/simulation.py
index b25b76f85..4b65b973a 100644
--- a/pybamm/simulation.py
+++ b/pybamm/simulation.py
@@ -255,7 +255,7 @@ class Simulation:
             parameterised_model = new_parameter_values.process_model(
                 new_model, inplace=False
             )
-            self.experiment_unique_steps_to_model[repr(op)] = parameterised_model
+            self.experiment_unique_steps_to_model[op.basic_repr()] = parameterised_model
 
         # Set up rest model if experiment has start times
         if self.experiment.initial_start_time:
@@ -778,8 +778,8 @@ class Simulation:
                     else:
                         dt = op_conds.duration
                     op_conds_str = str(op_conds)
-                    model = self.op_conds_to_built_models[repr(op_conds)]
-                    solver = self.op_conds_to_built_solvers[repr(op_conds)]
+                    model = self.op_conds_to_built_models[op_conds.basic_repr()]
+                    solver = self.op_conds_to_built_solvers[op_conds.basic_repr()]
 
                     logs["step number"] = (step_num, cycle_length)
                     logs["step operating conditions"] = op_conds_str
diff --git a/pybamm/solvers/solution.py b/pybamm/solvers/solution.py
index ecc8f7b70..4c9ccb993 100644
--- a/pybamm/solvers/solution.py
+++ b/pybamm/solvers/solution.py
@@ -443,7 +443,7 @@ class Solution(object):
 
     @initial_start_time.setter
     def initial_start_time(self, value):
-        """Updates the reason for termination"""
+        """Updates the initial start time of the experiment"""
         self._initial_start_time = value
 
     def set_summary_variables(self, all_summary_variables):
diff --git a/pybamm/step/_steps_util.py b/pybamm/step/_steps_util.py
index 879461b73..e524bc606 100644
--- a/pybamm/step/_steps_util.py
+++ b/pybamm/step/_steps_util.py
@@ -71,22 +71,25 @@ class _Step:
     ):
         self.type = typ
 
-        # Record all the args for repr
-        self.args = f"{typ}, {value}"
+        # Record all the args for repr and hash
+        self.repr_args = f"{typ}, {value}"
+        self.hash_args = f"{typ}, {value}"
         if duration:
-            self.args += f", duration={duration}"
+            self.repr_args += f", duration={duration}"
         if termination:
-            self.args += f", termination={termination}"
+            self.repr_args += f", termination={termination}"
+            self.hash_args += f", termination={termination}"
         if period:
-            self.args += f", period={period}"
+            self.repr_args += f", period={period}"
         if temperature:
-            self.args += f", temperature={temperature}"
+            self.repr_args += f", temperature={temperature}"
+            self.hash_args += f", temperature={temperature}"
         if tags:
-            self.args += f", tags={tags}"
+            self.repr_args += f", tags={tags}"
         if start_time:
-            self.args += f", start_time={start_time}"
+            self.repr_args += f", start_time={start_time}"
         if description:
-            self.args += f", description={description}"
+            self.repr_args += f", description={description}"
 
         # Check if drive cycle
         self.is_drive_cycle = isinstance(value, np.ndarray)
@@ -158,7 +161,15 @@ class _Step:
             return repr(self)
 
     def __repr__(self):
-        return f"_Step({self.args})"
+        return f"_Step({self.repr_args})"
+
+    def basic_repr(self):
+        """
+        Return a basic representation of the step, only with type, value, termination
+        and temperature, which are the variables involved in processing the model. Also
+        used for hashing.
+        """
+        return f"_Step({self.hash_args})"
 
     def to_dict(self):
         """
@@ -184,13 +195,11 @@ class _Step:
     def __eq__(self, other):
         return (
             isinstance(other, _Step)
-            and self.__repr__() == other.__repr__()
-            and self.next_start_time == other.next_start_time
-            and self.end_time == other.end_time
+            and self.hash_args == other.hash_args
         )
 
     def __hash__(self):
-        return hash(repr(self))
+        return hash(self.basic_repr())
 
     @property
     def unit(self):


Test patch:
diff --git a/tests/unit/test_experiments/test_simulation_with_experiment.py b/tests/unit/test_experiments/test_simulation_with_experiment.py
index 13531281d..6688fae5b 100644
--- a/tests/unit/test_experiments/test_simulation_with_experiment.py
+++ b/tests/unit/test_experiments/test_simulation_with_experiment.py
@@ -38,8 +38,12 @@ class TestSimulationExperiment(TestCase):
             [3600, 3 / Crate * 3600, 24 * 3600, 24 * 3600],
         )
 
-        model_I = sim.experiment_unique_steps_to_model[repr(op_conds[1])]  # CC charge
-        model_V = sim.experiment_unique_steps_to_model[repr(op_conds[2])]  # CV hold
+        model_I = sim.experiment_unique_steps_to_model[
+            op_conds[1].basic_repr()
+        ]  # CC charge
+        model_V = sim.experiment_unique_steps_to_model[
+            op_conds[2].basic_repr()
+        ]  # CV hold
         self.assertIn(
             "Current cut-off [A] [experiment]",
             [event.name for event in model_V.events],
@@ -729,6 +733,38 @@ class TestSimulationExperiment(TestCase):
         # test that the final time is correct (i.e. starting solution correctly set)
         self.assertEqual(new_solution["Time [s]"].entries[-1], 5400)
 
+    def test_experiment_start_time_identical_steps(self):
+        # Test that if we have the same step twice, with different start times,
+        # they get processed only once
+        model = pybamm.lithium_ion.SPM()
+
+        experiment = pybamm.Experiment(
+            [
+                pybamm.step.string(
+                    "Discharge at C/2 for 10 minutes",
+                    start_time=datetime(2023, 1, 1, 8, 0, 0),
+                ),
+                pybamm.step.string("Discharge at C/3 for 10 minutes"),
+                pybamm.step.string(
+                    "Discharge at C/2 for 10 minutes",
+                    start_time=datetime(2023, 1, 1, 10, 0, 0),
+                ),
+                pybamm.step.string("Discharge at C/3 for 10 minutes"),
+            ]
+        )
+
+        sim = pybamm.Simulation(model, experiment=experiment)
+        sim.solve(calc_esoh=False)
+
+        # Check that there are 4 steps
+        self.assertEqual(len(experiment.operating_conditions_steps), 4)
+
+        # Check that there are only 2 unique steps
+        self.assertEqual(len(sim.experiment.unique_steps), 2)
+
+        # Check that there are only 3 built models (unique steps + padding rest)
+        self.assertEqual(len(sim.op_conds_to_built_models), 3)
+
 
 if __name__ == "__main__":
     print("Add -v for more debug output")
diff --git a/tests/unit/test_parameters/test_bpx.py b/tests/unit/test_parameters/test_bpx.py
index 8288655b8..2559641d7 100644
--- a/tests/unit/test_parameters/test_bpx.py
+++ b/tests/unit/test_parameters/test_bpx.py
@@ -9,7 +9,6 @@ import json
 import pybamm
 import copy
 
-
 class TestBPX(TestCase):
     def setUp(self):
         self.base = {
@@ -161,7 +160,26 @@ class TestBPX(TestCase):
             json.dump(bpx_obj, tmp)
             tmp.flush()
 
-            pybamm.ParameterValues.create_from_bpx(tmp.name)
+            param = pybamm.ParameterValues.create_from_bpx(tmp.name)
+
+            # Function to check that functional parameters output constants
+            def check_constant_output(func):
+                stos = [0, 1]
+                T = 298.15
+                p_vals = [func(sto, T) for sto in stos]
+                self.assertEqual(p_vals[0], p_vals[1])
+
+            for electrode in ["Negative", "Positive"]:
+                D = param[f"{electrode} electrode diffusivity [m2.s-1]"]
+                dUdT = param[f"{electrode} electrode OCP entropic change [V.K-1]"]
+                check_constant_output(D)
+                check_constant_output(dUdT)
+
+            kappa = param["Electrolyte conductivity [S.m-1]"]
+            De = param["Electrolyte diffusivity [m2.s-1]"]
+            check_constant_output(kappa)
+            check_constant_output(De)
+
 
     def test_table_data(self):
         bpx_obj = copy.copy(self.base)
@@ -222,6 +240,66 @@ class TestBPX(TestCase):
         with self.assertRaisesRegex(ValueError, "Target SOC"):
             pybamm.ParameterValues.create_from_bpx("blah.json", target_soc=10)
 
+    def test_bpx_arrhenius(self):
+        bpx_obj = copy.copy(self.base)
+
+        filename = "tmp.json"
+        with tempfile.NamedTemporaryFile(
+            suffix=filename, delete=False, mode="w"
+        ) as tmp:
+            # write to a tempory file so we can
+            # get the source later on using inspect.getsource
+            # (as long as the file still exists)
+            json.dump(bpx_obj, tmp)
+            tmp.flush()
+
+            pv = pybamm.ParameterValues.create_from_bpx(tmp.name)
+
+
+        def arrhenius_assertion(pv, param_key, Ea_key):
+            sto = 0.5
+            T = 300
+            c_e = 1000
+            c_s_surf = 15000
+            c_s_max = 20000
+            T_ref = pv["Reference temperature [K]"]
+            Ea = pv[Ea_key]
+
+            if "exchange-current" in param_key:
+                eval_ratio = (
+                    pv[param_key](c_e, c_s_surf, c_s_max, T).value
+                    / pv[param_key](c_e, c_s_surf, c_s_max, T_ref).value
+                    )
+            else:
+                eval_ratio = (
+                    pv[param_key](sto, T).value
+                    / pv[param_key](sto, T_ref).value
+                )
+
+            calc_ratio = pybamm.exp(Ea / pybamm.constants.R * (1 / T_ref - 1 / T)).value
+
+            self.assertAlmostEqual(eval_ratio, calc_ratio)
+
+        param_keys = [
+            "Electrolyte conductivity [S.m-1]",
+            "Electrolyte diffusivity [m2.s-1]",
+            "Negative electrode diffusivity [m2.s-1]",
+            "Positive electrode diffusivity [m2.s-1]",
+            "Positive electrode exchange-current density [A.m-2]",
+            "Negative electrode exchange-current density [A.m-2]",
+        ]
+
+        Ea_keys = [
+            "Electrolyte conductivity activation energy [J.mol-1]",
+            "Electrolyte diffusivity activation energy [J.mol-1]",
+            "Negative electrode diffusivity activation energy [J.mol-1]",
+            "Positive electrode diffusivity activation energy [J.mol-1]",
+            "Positive electrode reaction rate constant activation energy [J.mol-1]",
+            "Negative electrode reaction rate constant activation energy [J.mol-1]",
+        ]
+
+        for param_key, Ea_key in zip(param_keys, Ea_keys):
+            arrhenius_assertion(pv, param_key, Ea_key)
 
 if __name__ == "__main__":
     print("Add -v for more debug output")