diff --git a/demos/helheim/basal_friction.py b/demos/helheim/basal_friction.py
index 91617d5..515bff0 100644
--- a/demos/helheim/basal_friction.py
+++ b/demos/helheim/basal_friction.py
@@ -156,7 +156,7 @@ def regularization(u):
 
 def bed_friction(**kwargs):
     u, q = map(kwargs.get, ("velocity", "log_friction"))
-    C = τ_avg / u_avg ** (1 / m) * exp(-q)
+    C = Constant(τ_avg) / Constant(u_avg) ** (1 / m) * exp(-q)
     return icepack.models.friction.bed_friction(velocity=u, friction=C)
 
 # Compute an initial estimate for the ice velocity
@@ -194,7 +194,7 @@ def regularization(q):
     α = Constant(args.regularization)
     return 0.5 * α**2 / Ω * inner(grad(q), grad(q)) * dx
 
-# Set up an inverse solver / problem
+# Solve a statistical estimation problem for the log-friction
 problem = icepack.statistics.StatisticsProblem(
     simulation=simulation,
     loss_functional=loss_functional,
@@ -211,9 +211,9 @@ def regularization(q):
 
 q = estimator.solve()
 
-import matplotlib.pyplot as plt
-fig, axes = plt.subplots()
-axes.set_aspect("equal")
-colors = firedrake.tripcolor(q, axes=axes)
-fig.colorbar(colors)
-plt.show()
+# Save the result to disk
+with firedrake.CheckpointFile(args.output, "w") as chk:
+    chk.save_function(q, name="log_friction")
+    chk.save_function(u, name="velocity")
+    chk.h5pyfile.attrs["mean_stress"] = τ_avg
+    chk.h5pyfile.attrs["mean_speed"] = u_avg