PQN fixes from David's fix_pqn branch; added strong wolfe line search…

… to PQN for BFGS equivalency;cleaned up complex tests from SPG;complex function test now passes with PQN

math/src/main/scala/breeze/optimize/ProjectedQuasiNewton.scala

@@ -62,8 +62,6 @@ class CompactHessian(M: DenseMatrix[Double], Y: RingBuffer[DenseVector[Double]],
   lazy val N = DenseMatrix.horzcat(collectionOfVectorsToMatrix(S).t * sigma, collectionOfVectorsToMatrix(Y).t)
 }
 
-//TO DO:
-//1. Modify SpectralProjectedGradient to implement SpaRSA so that it can handle proximal algorithms
 class ProjectedQuasiNewton(tolerance: Double = 1e-6,
                            val m: Int = 10,
                            val initFeas: Boolean = false,
@@ -109,45 +107,27 @@ class ProjectedQuasiNewton(tolerance: Double = 1e-6,
     }
   }
 
-  protected def determineStepSize(state: State, fn: DiffFunction[DenseVector[Double]], dir: DenseVector[Double]): Double = {
-    if (state.iter == 0)
-      return scala.math.min(1.0, 1.0 / norm(state.grad,1.0))
-    val dirnorm = norm(dir, Double.PositiveInfinity)
-    if(dirnorm < 1E-10) return 0.0
-    import state._
-    // Backtracking line-search
-    var accepted = false
-    var lambda = 1.0
-    val gTd = grad dot dir
-    var srchit = 0
-
-    do {
-      val candx = x + dir * lambda
-      val candf = fn.valueAt(candx)
-      val suffdec = gamma * lambda * gTd
-
-      if (testOpt && srchit > 0) {
-        logger.debug(f"PQN:    SrchIt $srchit%4d: f $candf%-10.4f t $lambda%-10.4f\n")
-      }
-
-      if (candf < state.adjustedValue + suffdec) {
-        accepted = true
-      } else if (srchit >= maxSrchIt) {
-        accepted = true
-      } else {
-        lambda *= 0.5
-        srchit = srchit + 1
-      }
-    } while (!accepted)
-
-    if (srchit >= maxSrchIt) {
-      logger.info("PQN: Line search cannot make further progress")
-      throw new LineSearchFailed(norm(state.grad,Double.PositiveInfinity), norm(dir, Double.PositiveInfinity))
-    }
-    lambda
+  /**
+   * Given a direction, perform a Strong Wolfe Line Search
+   *
+   * @param state the current state
+   * @param f The objective
+   * @param dir The step direction
+   * @return stepSize
+   */
+  protected def determineStepSize(state: State, f: DiffFunction[DenseVector[Double]], dir: DenseVector[Double]) = {
+    val x = state.x
+    val grad = state.grad
+
+    val ff = LineSearch.functionFromSearchDirection(f, x, dir)
+    val search = new StrongWolfeLineSearch(maxZoomIter = 10, maxLineSearchIter = maxSrchIt) // TODO: Need good default values here.
+    val alpha = search.minimize(ff, if(state.iter == 0.0) min(1.0, 1.0/norm(dir)) else 1.0)
+
+    if(alpha * norm(grad) < 1E-10)
+      throw new StepSizeUnderflow
+    alpha
   }
 
-
   protected def takeStep(state: State, dir: DenseVector[Double], stepSize: Double): DenseVector[Double] = {
     projection(state.x + dir * stepSize)
   }

math/src/test/scala/breeze/optimize/SpectralProjectedGradientTest.scala

@@ -76,42 +76,6 @@ class SpectralProjectedGradientTest extends PropSpec with PropertyChecks with Op
     }
   }
 
-
-  property("optimize a simple multivariate gaussian with l2 regularization") {
-    val optimizer = new SpectralProjectedGradient[DenseVector[Double]](tolerance = 1.0E-5)
-
-    forAll { init: DenseVector[Double] =>
-      val f = new DiffFunction[DenseVector[Double]] {
-        def calculate(x: DenseVector[Double]) = {
-          (norm((x - 3.0) :^ 2.0, 1), (x * 2.0) - 6.0)
-        }
-      }
-
-      val targetValue = 3 / (1.0 / 2 + 1)
-      val result = optimizer.minimize(DiffFunction.withL2Regularization(f, 1.0), init)
-      result should beSimilarTo(DenseVector.ones[Double](init.size) * targetValue, allowedDeviation = 3E-3 * result.size)
-    }
-  }
-
-  property("optimize a complicated function without projection") {
-    val optimizer = new SpectralProjectedGradient[DenseVector[Double]](tolerance = 1.0E-5)
-
-    forAll { a: DenseVector[Double] =>
-      whenever(min(a) >= -3.0 && max(a) <= 3.0) {
-        val init = DenseVector.rand(a.size)
-        val f = new DiffFunction[DenseVector[Double]] {
-          def calculate(x: DenseVector[Double]) = {
-            (sum(exp((x :^ 2.0) :- (a :* x))), (x * 2.0 :- a) :* exp(x :^ 2.0 :- a :* x))
-          }
-        }
-
-        val result = optimizer.minimize(f, init)
-        val minimum = f(a / 2.0)
-        f(result) should be(minimum +- abs(minimum) * 1E-2)
-      }
-    }
-  }
-
   property("simple linear solve without projection") {
     val n = 5
     val H = new DenseMatrix(n, n, Array(1.8984250861699135,0.5955576666769438,-1.484430453342902,-1.0434994471390804,-3.675310432634351,0.5955576666769438,0.9090751938470876,-2.146380947361661,-0.13037609428980368,-0.40639564652095117,-1.484430453342902,-2.146380947361661,10.262733520770384,-6.097698907163584,2.29625304115155,-1.0434994471390804,-0.13037609428980368,-6.097698907163584,27.775920405610677,-5.574220233644466,-3.675310432634351,-0.40639564652095117,2.29625304115155,-5.574220233644466,12.21329172136971))