Use more references

src/samplefunctional.hh

@@ -39,7 +39,8 @@ template <int dim> class SampleFunctional {
     // negative_projection() divides by it
     if (x.two_norm2() == 0.0) {
       // If there is a direction of descent, this is it
-      SmallVector const d = smoothGradient(x);
+      SmallVector d;
+      smoothGradient(x, d);
 
       Interval<double> D;
       phi.directionalSubDiff(x, d, D);
@@ -56,8 +57,10 @@ template <int dim> class SampleFunctional {
       return;
     }
 
-    SmallVector const pg = upperGradient(x);
-    SmallVector const mg = lowerGradient(x);
+    SmallVector pg;
+    upperGradient(x, pg);
+    SmallVector mg;
+    lowerGradient(x, mg);
     double const pgx = pg * x;
     double const mgx = mg * x;
     if (pgx >= 0 && mgx >= 0) {
@@ -73,7 +76,9 @@ template <int dim> class SampleFunctional {
     } else {
       // Includes the case that pg points in direction x and mg
       // points in direction -x. The projection will then be zero.
-      ret = negative_projection(smoothGradient(x), x);
+      SmallVector d;
+      smoothGradient(x, d);
+      ret = negative_projection(d, x);
       dverb << "## Directional derivative (as per scalar product w/ "
                "semigradient): " << -(ret * ret)
             << " (coordinates of the restriction)" << std::endl;
@@ -87,25 +92,23 @@ template <int dim> class SampleFunctional {
 
 private:
   // Gradient of the smooth part
-  SmallVector smoothGradient(SmallVector const x) const {
-    SmallVector y;
+  void smoothGradient(SmallVector const &x, SmallVector &y) const {
     A.mv(x, y); // y = Av
     y -= b;     // y = Av - b
-    return y;
   }
 
-  SmallVector upperGradient(SmallVector const x) const {
-    SmallVector y;
+  void upperGradient(SmallVector const &x, SmallVector &y) const {
     phi.upperGradient(x, y);
-    y += smoothGradient(x);
-    return y;
+    SmallVector z;
+    smoothGradient(x, z);
+    y += z;
   }
 
-  SmallVector lowerGradient(SmallVector const x) const {
-    SmallVector y;
+  void lowerGradient(SmallVector const &x, SmallVector &y) const {
     phi.lowerGradient(x, y);
-    y += smoothGradient(x);
-    return y;
+    SmallVector z;
+    smoothGradient(x, z);
+    y += z;
   }
 
   // No normalising is done!