diff --git a/src/darsia/measure/wasserstein.py b/src/darsia/measure/wasserstein.py
index 73102614..05ebdf52 100644
--- a/src/darsia/measure/wasserstein.py
+++ b/src/darsia/measure/wasserstein.py
@@ -705,8 +705,8 @@ def vector_face_flux_norm(self, flat_flux: np.ndarray, mode: str) -> np.ndarray:
                             subcell_flux, 2, axis=-1
                         ).ravel("F")[cells]
 
-                # Average over the subcells using harmonic averaging
-                flat_flux_norm[faces] = hmean(subcell_flux_norm[faces], axis=1)
+            # Average over the subcells using harmonic averaging
+            flat_flux_norm = hmean(subcell_flux_norm, axis=1)
 
         elif mode == "face_based":
             if not hasattr(self, "face_reconstruction"):