diff --git a/examples/scripts/compare_lead_acid.py b/examples/scripts/compare_lead_acid.py
index 29824fada9..91def5ce74 100644
--- a/examples/scripts/compare_lead_acid.py
+++ b/examples/scripts/compare_lead_acid.py
@@ -18,15 +18,15 @@
 
 # load models
 models = [
-    # pybamm.lead_acid.LOQS(
-    #     {"surface form": "differential", "bc_options": {"dimensionality": 1}},
-    #     name="3D LOQS model",
-    # ),
-    pybamm.lead_acid.LOQS(),
+    pybamm.lead_acid.LOQS(
+        {"surface form": "differential", "bc_options": {"dimensionality": 1}},
+        name="3D LOQS model",
+    ),
+    # pybamm.lead_acid.LOQS(),
     # pybamm.lead_acid.FOQS(),
-    pybamm.lead_acid.CompositeExtended(),
+    # pybamm.lead_acid.CompositeExtended(),
     # # pybamm.lead_acid.Composite({"surface form": "algebraic"}),
-    pybamm.lead_acid.NewmanTiedemann(),
+    # pybamm.lead_acid.NewmanTiedemann(),
 ]
 
 # load parameter values and process models and geometry
@@ -59,22 +59,22 @@
 
 # solve model
 solutions = [None] * len(models)
-t_eval = np.linspace(0, 5, 100)
+t_eval = np.linspace(0, 1, 100)
 for i, model in enumerate(models):
     solution = model.default_solver.solve(model, t_eval)
     solutions[i] = solution
 
 # plot
 output_variables = [
-    [
-        "Average negative electrode interfacial current density [A.m-2]",
-        "Average positive electrode interfacial current density [A.m-2]",
-    ],
-    "Average negative electrode surface potential difference [V]",
-    "Average positive electrode surface potential difference [V]",
-    "Electrolyte concentration",
-    "Electrolyte flux",
-    "Terminal voltage [V]",
+    # [
+    #     "Average negative electrode interfacial current density [A.m-2]",
+    #     "Average positive electrode interfacial current density [A.m-2]",
+    # ],
+    # "Average negative electrode surface potential difference [V]",
+    # "Average positive electrode surface potential difference [V]",
+    # "Electrolyte concentration",
+    # "Electrolyte flux",
+    "Terminal voltage [V]"
 ]
 plot = pybamm.QuickPlot(models, mesh, solutions, output_variables)
 plot.dynamic_plot()
diff --git a/pybamm/models/full_battery_models/base_battery_model.py b/pybamm/models/full_battery_models/base_battery_model.py
index 1981f886ca..b42c6629f4 100644
--- a/pybamm/models/full_battery_models/base_battery_model.py
+++ b/pybamm/models/full_battery_models/base_battery_model.py
@@ -433,11 +433,17 @@ def set_voltage_variables(self):
             V = pybamm.BoundaryValue(phi_s_p, "right")
             V_dim = pybamm.BoundaryValue(phi_s_p_dim, "right")
         elif self.options["bc_options"]["dimensionality"] == 1:
-            # TO DO: add terminal voltage in 1plus1D
-            phi_s_p = self.variables["Positive electrode potential"]
-            phi_s_p_dim = self.variables["Positive electrode potential [V]"]
-            V = pybamm.BoundaryValue(phi_s_p, "right")
-            V_dim = pybamm.BoundaryValue(phi_s_p_dim, "right")
+            phi_s_cn = self.variables["Negative current collector potential"]
+            phi_s_cp = self.variables["Positive current collector potential"]
+            phi_s_cn_dim = self.variables["Negative current collector potential [V]"]
+            phi_s_cp_dim = self.variables["Positive current collector potential [V]"]
+            # In 1D both tabs are at "right"
+            V = pybamm.BoundaryValue(phi_s_cp, "right") - pybamm.BoundaryValue(
+                phi_s_cn, "right"
+            )
+            V_dim = pybamm.BoundaryValue(phi_s_cp_dim, "right") - pybamm.BoundaryValue(
+                phi_s_cn_dim, "right"
+            )
         elif self.options["bc_options"]["dimensionality"] == 2:
             phi_s_cn = self.variables["Negative current collector potential"]
             phi_s_cp = self.variables["Positive current collector potential"]