diff --git a/src/localgeodesicfestiffness.hh b/src/localgeodesicfestiffness.hh
index 9ffddea3c7ddf3a96e0097aa2beec906a4f8f5d6..d21b23cb2357a300e213f3b8f9635a9ae208d858 100644
--- a/src/localgeodesicfestiffness.hh
+++ b/src/localgeodesicfestiffness.hh
@@ -396,9 +396,11 @@ assembleGradient(const Entity& element,
         
         for (int j=0; j<blocksize; j++) {
             
-            // Brute force: the return value does not have unit norm.  Stuff it in there anyways
-            forwardSolution[i].data_  = infinitesimalVariation(localSolution[i],  eps, j);
-            backwardSolution[i].data_ = infinitesimalVariation(localSolution[i], -eps, j);
+            // The return value does not have unit norm.  But assigning it to a UnitVector object
+            // will normalize it.  This amounts to an extension of the energy functional 
+            // to a neighborhood around S^n
+            forwardSolution[i]  = infinitesimalVariation(localSolution[i],  eps, j);
+            backwardSolution[i] = infinitesimalVariation(localSolution[i], -eps, j);
 
             localGradient[i][j] = (energy(element,forwardSolution) - energy(element,backwardSolution))
                 / (2*eps);