commons-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From er...@apache.org
Subject svn commit: r1188383 [2/3] - in /commons/proper/math/trunk/src: main/java/org/apache/commons/math/optimization/direct/BOBYQAOptimizer.java test/java/org/apache/commons/math/optimization/direct/BOBYQAOptimizerTest.java
Date Mon, 24 Oct 2011 21:09:32 GMT

Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/direct/BOBYQAOptimizer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/direct/BOBYQAOptimizer.java?rev=1188383&r1=1188382&r2=1188383&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/direct/BOBYQAOptimizer.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/optimization/direct/BOBYQAOptimizer.java Mon Oct 24 21:09:32 2011
@@ -113,6 +113,66 @@ public class BOBYQAOptimizer
      * Index of the interpolation point at the trust region center.
      */
     private int trustRegionCenterInterpolationPointIndex;
+    /**
+     * XXX "bmat" in the original code.
+     */
+    private Array2DRowRealMatrix bMatrix;
+    /**
+     * XXX "zmat" in the original code.
+     */
+    private Array2DRowRealMatrix zMatrix;
+    /**
+     * XXX "xpt" in the original code.
+     */
+    private Array2DRowRealMatrix interpolationPoints;
+    /**
+     * XXX "xbase" in the original code.
+     */
+    private ArrayRealVector originShift;
+    /**
+     * XXX "fval" in the original code.
+     */
+    private ArrayRealVector fAtInterpolationPoints;
+    /**
+     * XXX "xopt" in the original code.
+     */
+    private ArrayRealVector trustRegionCenterOffset;
+    /**
+     * XXX "gopt" in the original code.
+     */
+    private ArrayRealVector gradientAtTrustRegionCenter;
+    /**
+     * XXX "sl" in the original code.
+     */
+    private ArrayRealVector lowerDifference;
+    /**
+     * XXX "su" in the original code.
+     */
+    private ArrayRealVector upperDifference;
+    /**
+     * XXX "pq" in the original code.
+     */
+    private ArrayRealVector modelSecondDerivativesParameters;
+    /**
+     * XXX "xnew" in the original code.
+     */
+    private ArrayRealVector newPoint;
+    /**
+     * XXX "xalt" in the original code.
+     */
+    private ArrayRealVector alternativeNewPoint;
+    /**
+     * XXX "d__" in the original code.
+     */
+    private ArrayRealVector trialStepPoint;
+    /**
+     * XXX "vlag" in the original code.
+     */
+    private ArrayRealVector lagrangeValuesAtNewPoint;
+    /**
+     * XXX "hq" in the original code.
+     */
+    private ArrayRealVector modelSecondDerivativesValues;
 
     /**
      * @param numberOfInterpolationPoints Number of interpolation conditions.
@@ -212,34 +272,9 @@ public class BOBYQAOptimizer
      * @return
      */
     private double bobyqa() {
-        // System.out.println("bobyqa"); // XXX
+        printMethod(); // XXX
 
         final int n = currentBest.getDimension();
-        final int npt = numberOfInterpolationPoints;
-
-        final int np = n + 1;
-        final int ndim = npt + n;
-
-        // Partition the working space array, so that different parts of it can
-        // be treated separately during the calculation of BOBYQB. The partition
-        // requires the first (NPT+2)*(NPT+N)+3*N*(N+5)/2 elements of W plus the
-        // space that is taken by the last array in the argument list of BOBYQB.
-
-        final ArrayRealVector xbase = new ArrayRealVector(n);
-        final Array2DRowRealMatrix xpt = new Array2DRowRealMatrix(npt, n);
-        final ArrayRealVector fval = new ArrayRealVector(npt);
-        final ArrayRealVector xopt = new ArrayRealVector(n);
-        final ArrayRealVector gopt = new ArrayRealVector(n);
-        final ArrayRealVector hq = new ArrayRealVector(n * np / 2);
-        final ArrayRealVector pq = new ArrayRealVector(npt);
-        final Array2DRowRealMatrix bmat = new Array2DRowRealMatrix(ndim, n);
-        final Array2DRowRealMatrix zmat = new Array2DRowRealMatrix(npt, (npt - np));
-        final ArrayRealVector sl = new ArrayRealVector(n);
-        final ArrayRealVector su = new ArrayRealVector(n);
-        final ArrayRealVector xnew = new ArrayRealVector(n);
-        final ArrayRealVector xalt = new ArrayRealVector(n);
-        final ArrayRealVector d__ = new ArrayRealVector(n);
-        final ArrayRealVector vlag = new ArrayRealVector(ndim);
 
         // Return if there is insufficient space between the bounds. Modify the
         // initial X if necessary in order to avoid conflicts between the bounds
@@ -250,53 +285,39 @@ public class BOBYQAOptimizer
 
         for (int j = 0; j < n; j++) {
             final double boundDiff = boundDifference[j];
-            sl.setEntry(j, lowerBound[j] - currentBest.getEntry(j));
-            su.setEntry(j, upperBound[j] - currentBest.getEntry(j));
-            if (sl.getEntry(j) >= -initialTrustRegionRadius) {
-                if (sl.getEntry(j) >= ZERO) {
+            lowerDifference.setEntry(j, lowerBound[j] - currentBest.getEntry(j));
+            upperDifference.setEntry(j, upperBound[j] - currentBest.getEntry(j));
+            if (lowerDifference.getEntry(j) >= -initialTrustRegionRadius) {
+                if (lowerDifference.getEntry(j) >= ZERO) {
                     currentBest.setEntry(j, lowerBound[j]);
-                    sl.setEntry(j, ZERO);
-                    su.setEntry(j, boundDiff);
+                    lowerDifference.setEntry(j, ZERO);
+                    upperDifference.setEntry(j, boundDiff);
                 } else {
                     currentBest.setEntry(j, lowerBound[j] + initialTrustRegionRadius);
-                    sl.setEntry(j, -initialTrustRegionRadius);
+                    lowerDifference.setEntry(j, -initialTrustRegionRadius);
                     // Computing MAX
                     final double deltaOne = upperBound[j] - currentBest.getEntry(j);
-                    su.setEntry(j, Math.max(deltaOne, initialTrustRegionRadius));
+                    upperDifference.setEntry(j, Math.max(deltaOne, initialTrustRegionRadius));
                 }
-            } else if (su.getEntry(j) <= initialTrustRegionRadius) {
-                if (su.getEntry(j) <= ZERO) {
+            } else if (upperDifference.getEntry(j) <= initialTrustRegionRadius) {
+                if (upperDifference.getEntry(j) <= ZERO) {
                     currentBest.setEntry(j, upperBound[j]);
-                    sl.setEntry(j, -boundDiff);
-                    su.setEntry(j, ZERO);
+                    lowerDifference.setEntry(j, -boundDiff);
+                    upperDifference.setEntry(j, ZERO);
                 } else {
                     currentBest.setEntry(j, upperBound[j] - initialTrustRegionRadius);
                     // Computing MIN
                     final double deltaOne = lowerBound[j] - currentBest.getEntry(j);
                     final double deltaTwo = -initialTrustRegionRadius;
-                    sl.setEntry(j, Math.min(deltaOne, deltaTwo));
-                    su.setEntry(j, initialTrustRegionRadius);
+                    lowerDifference.setEntry(j, Math.min(deltaOne, deltaTwo));
+                    upperDifference.setEntry(j, initialTrustRegionRadius);
                 }
             }
         }
 
         // Make the call of BOBYQB.
 
-        return bobyqb(xbase,
-                      xpt,
-                      fval,
-                      xopt,
-                      gopt,
-                      hq,
-                      pq,
-                      bmat,
-                      zmat,
-                      sl,
-                      su,
-                      xnew,
-                      xalt,
-                      d__,
-                      vlag);
+        return bobyqb();
     } // bobyqa
 
     // ----------------------------------------------------------------------------------------
@@ -334,41 +355,10 @@ public class BOBYQAOptimizer
      *     W is a one-dimensional array that is used for working space. Its length
      *       must be at least 3*NDIM = 3*(NPT+N).
      *
-     * @param xbase
-     * @param xpt
-     * @param fval
-     * @param xopt
-     * @param gopt
-     * @param hq
-     * @param pq
-     * @param bmat
-     * @param zmat
-     * @param sl
-     * @param su
-     * @param xnew
-     * @param xalt
-     * @param d__
-     * @param vlag
      * @return
      */
-    private double bobyqb(
-            ArrayRealVector xbase,
-            Array2DRowRealMatrix xpt,
-            ArrayRealVector fval,
-            ArrayRealVector xopt,
-            ArrayRealVector gopt,
-            ArrayRealVector hq,
-            ArrayRealVector pq,
-            Array2DRowRealMatrix bmat,
-            Array2DRowRealMatrix zmat,
-            ArrayRealVector sl,
-            ArrayRealVector su,
-            ArrayRealVector xnew,
-            ArrayRealVector xalt,
-            ArrayRealVector d__,
-            ArrayRealVector vlag
-    ) {
-        // System.out.println("bobyqb"); // XXX
+    private double bobyqb() {
+        printMethod(); // XXX
 
         final int n = currentBest.getDimension();
         final int npt = numberOfInterpolationPoints;
@@ -399,17 +389,15 @@ public class BOBYQAOptimizer
 
         trustRegionCenterInterpolationPointIndex = 0;
 
-        prelim(currentBest, xbase,
-               xpt, fval, gopt, hq, pq, bmat,
-               zmat, sl, su);
+        prelim();
         double xoptsq = ZERO;
         for (int i = 0; i < n; i++) {
-            xopt.setEntry(i, xpt.getEntry(trustRegionCenterInterpolationPointIndex, i));
+            trustRegionCenterOffset.setEntry(i, interpolationPoints.getEntry(trustRegionCenterInterpolationPointIndex, i));
             // Computing 2nd power
-            final double deltaOne = xopt.getEntry(i);
+            final double deltaOne = trustRegionCenterOffset.getEntry(i);
             xoptsq += deltaOne * deltaOne;
         }
-        double fsave = fval.getEntry(0);
+        double fsave = fAtInterpolationPoints.getEntry(0);
         final int kbase = 0;
 
         // Complete the settings that are required for the iterative procedure.
@@ -439,14 +427,15 @@ public class BOBYQAOptimizer
         int state = 20;
         for(;;) switch (state) {
         case 20: {
+            printState(20); // XXX
             if (trustRegionCenterInterpolationPointIndex != kbase) {
                 int ih = 0;
                 for (int j = 0; j < n; j++) {
                     for (int i = 0; i <= j; i++) {
                         if (i < j) {
-                            gopt.setEntry(j,  gopt.getEntry(j) + hq.getEntry(ih) * xopt.getEntry(i));
+                            gradientAtTrustRegionCenter.setEntry(j,  gradientAtTrustRegionCenter.getEntry(j) + modelSecondDerivativesValues.getEntry(ih) * trustRegionCenterOffset.getEntry(i));
                         }
-                        gopt.setEntry(i,  gopt.getEntry(i) + hq.getEntry(ih) * xopt.getEntry(j));
+                        gradientAtTrustRegionCenter.setEntry(i,  gradientAtTrustRegionCenter.getEntry(i) + modelSecondDerivativesValues.getEntry(ih) * trustRegionCenterOffset.getEntry(j));
                         ih++;
                     }
                 }
@@ -454,11 +443,11 @@ public class BOBYQAOptimizer
                     for (int k = 0; k < npt; k++) {
                         double temp = ZERO;
                         for (int j = 0; j < n; j++) {
-                            temp += xpt.getEntry(k, j) * xopt.getEntry(j);
+                            temp += interpolationPoints.getEntry(k, j) * trustRegionCenterOffset.getEntry(j);
                         }
-                        temp *= pq.getEntry(k);
+                        temp *= modelSecondDerivativesParameters.getEntry(k);
                         for (int i = 0; i < n; i++) {
-                            gopt.setEntry(i, gopt.getEntry(i) + temp * xpt.getEntry(k, i));
+                            gradientAtTrustRegionCenter.setEntry(i, gradientAtTrustRegionCenter.getEntry(i) + temp * interpolationPoints.getEntry(k, i));
                         }
                     }
                     throw new PathIsExploredException(); // XXX
@@ -474,14 +463,14 @@ public class BOBYQAOptimizer
 
         }
         case 60: {
+            printState(60); // XXX
             final ArrayRealVector gnew = new ArrayRealVector(n);
             final ArrayRealVector xbdi = new ArrayRealVector(n);
             final ArrayRealVector s = new ArrayRealVector(n);
             final ArrayRealVector hs = new ArrayRealVector(n);
             final ArrayRealVector hred = new ArrayRealVector(n);
 
-            final double[] dsqCrvmin = trsbox(xpt, xopt, gopt, hq, pq, sl,
-                                              su, delta, xnew, d__, gnew, xbdi, s,
+            final double[] dsqCrvmin = trsbox(delta, gnew, xbdi, s,
                                               hs, hred);
             dsq = dsqCrvmin[0];
             crvmin = dsqCrvmin[1];
@@ -516,18 +505,18 @@ public class BOBYQAOptimizer
                 final double bdtol = errbig / rho;
                 for (int j = 0; j < n; j++) {
                     double bdtest = bdtol;
-                    if (xnew.getEntry(j) == sl.getEntry(j)) {
+                    if (newPoint.getEntry(j) == lowerDifference.getEntry(j)) {
                         bdtest = work1.getEntry(j);
                     }
-                    if (xnew.getEntry(j) == su.getEntry(j)) {
+                    if (newPoint.getEntry(j) == upperDifference.getEntry(j)) {
                         bdtest = -work1.getEntry(j);
                     }
                     if (bdtest < bdtol) {
-                        double curv = hq.getEntry((j + j * j) / 2);
+                        double curv = modelSecondDerivativesValues.getEntry((j + j * j) / 2);
                         for (int k = 0; k < npt; k++) {
                             // Computing 2nd power
-                            final double d1 = xpt.getEntry(k, j);
-                            curv += pq.getEntry(k) * (d1 * d1);
+                            final double d1 = interpolationPoints.getEntry(k, j);
+                            curv += modelSecondDerivativesParameters.getEntry(k) * (d1 * d1);
                         }
                         bdtest += HALF * curv * rho;
                         if (bdtest < bdtol) {
@@ -548,29 +537,30 @@ public class BOBYQAOptimizer
 
         }
         case 90: {
+            printState(90); // XXX
             if (dsq <= xoptsq * ONE_OVER_A_THOUSAND) {
                 final double fracsq = xoptsq * ONE_OVER_FOUR;
                 double sumpq = ZERO;
                 // final RealVector sumVector
-                //     = new ArrayRealVector(npt, -HALF * xoptsq).add(xpt.operate(xopt));
+                //     = new ArrayRealVector(npt, -HALF * xoptsq).add(interpolationPoints.operate(trustRegionCenter));
                 for (int k = 0; k < npt; k++) {
-                    sumpq += pq.getEntry(k);
+                    sumpq += modelSecondDerivativesParameters.getEntry(k);
                     double sum = -HALF * xoptsq;
                     for (int i = 0; i < n; i++) {
-                        sum += xpt.getEntry(k, i) * xopt.getEntry(i);
+                        sum += interpolationPoints.getEntry(k, i) * trustRegionCenterOffset.getEntry(i);
                     }
                     // sum = sumVector.getEntry(k); // XXX "testAckley" and "testDiffPow" fail.
                     work2.setEntry(k, sum);
                     final double temp = fracsq - HALF * sum;
                     for (int i = 0; i < n; i++) {
-                        work1.setEntry(i, bmat.getEntry(k, i));
-                        vlag.setEntry(i, sum * xpt.getEntry(k, i) + temp * xopt.getEntry(i));
+                        work1.setEntry(i, bMatrix.getEntry(k, i));
+                        lagrangeValuesAtNewPoint.setEntry(i, sum * interpolationPoints.getEntry(k, i) + temp * trustRegionCenterOffset.getEntry(i));
                         final int ip = npt + i;
                         for (int j = 0; j <= i; j++) {
-                            bmat.setEntry(ip, j,
-                                          bmat.getEntry(ip, j)
-                                          + work1.getEntry(i) * vlag.getEntry(j)
-                                          + vlag.getEntry(i) * work1.getEntry(j));
+                            bMatrix.setEntry(ip, j,
+                                          bMatrix.getEntry(ip, j)
+                                          + work1.getEntry(i) * lagrangeValuesAtNewPoint.getEntry(j)
+                                          + lagrangeValuesAtNewPoint.getEntry(i) * work1.getEntry(j));
                         }
                     }
                 }
@@ -581,28 +571,28 @@ public class BOBYQAOptimizer
                     double sumz = ZERO;
                     double sumw = ZERO;
                     for (int k = 0; k < npt; k++) {
-                        sumz += zmat.getEntry(k, m);
-                        vlag.setEntry(k, work2.getEntry(k) * zmat.getEntry(k, m));
-                        sumw += vlag.getEntry(k);
+                        sumz += zMatrix.getEntry(k, m);
+                        lagrangeValuesAtNewPoint.setEntry(k, work2.getEntry(k) * zMatrix.getEntry(k, m));
+                        sumw += lagrangeValuesAtNewPoint.getEntry(k);
                     }
                     for (int j = 0; j < n; j++) {
-                        double sum = (fracsq * sumz - HALF * sumw) * xopt.getEntry(j);
+                        double sum = (fracsq * sumz - HALF * sumw) * trustRegionCenterOffset.getEntry(j);
                         for (int k = 0; k < npt; k++) {
-                            sum += vlag.getEntry(k) * xpt.getEntry(k, j);
+                            sum += lagrangeValuesAtNewPoint.getEntry(k) * interpolationPoints.getEntry(k, j);
                         }
                         work1.setEntry(j, sum);
                         for (int k = 0; k < npt; k++) {
-                            bmat.setEntry(k, j,
-                                          bmat.getEntry(k, j)
-                                          + sum * zmat.getEntry(k, m));
+                            bMatrix.setEntry(k, j,
+                                          bMatrix.getEntry(k, j)
+                                          + sum * zMatrix.getEntry(k, m));
                         }
                     }
                     for (int i = 0; i < n; i++) {
                         final int ip = i + npt;
                         final double temp = work1.getEntry(i);
                         for (int j = 0; j <= i; j++) {
-                            bmat.setEntry(ip, j,
-                                          bmat.getEntry(ip, j)
+                            bMatrix.setEntry(ip, j,
+                                          bMatrix.getEntry(ip, j)
                                           + temp * work1.getEntry(j));
                         }
                     }
@@ -613,26 +603,26 @@ public class BOBYQAOptimizer
 
                 int ih = 0;
                 for (int j = 0; j < n; j++) {
-                    work1.setEntry(j, -HALF * sumpq * xopt.getEntry(j));
+                    work1.setEntry(j, -HALF * sumpq * trustRegionCenterOffset.getEntry(j));
                     for (int k = 0; k < npt; k++) {
-                        work1.setEntry(j, work1.getEntry(j) + pq.getEntry(k) * xpt.getEntry(k, j));
-                        xpt.setEntry(k, j, xpt.getEntry(k, j) - xopt.getEntry(j));
+                        work1.setEntry(j, work1.getEntry(j) + modelSecondDerivativesParameters.getEntry(k) * interpolationPoints.getEntry(k, j));
+                        interpolationPoints.setEntry(k, j, interpolationPoints.getEntry(k, j) - trustRegionCenterOffset.getEntry(j));
                     }
                     for (int i = 0; i <= j; i++) {
-                         hq.setEntry(ih,
-                                    hq.getEntry(ih)
-                                    + work1.getEntry(i) * xopt.getEntry(j)
-                                    + xopt.getEntry(i) * work1.getEntry(j));
-                        bmat.setEntry(npt + i, j, bmat.getEntry(npt + j, i));
+                         modelSecondDerivativesValues.setEntry(ih,
+                                    modelSecondDerivativesValues.getEntry(ih)
+                                    + work1.getEntry(i) * trustRegionCenterOffset.getEntry(j)
+                                    + trustRegionCenterOffset.getEntry(i) * work1.getEntry(j));
+                        bMatrix.setEntry(npt + i, j, bMatrix.getEntry(npt + j, i));
                         ih++;
                     }
                 }
                 for (int i = 0; i < n; i++) {
-                    xbase.setEntry(i, xbase.getEntry(i) + xopt.getEntry(i));
-                    xnew.setEntry(i, xnew.getEntry(i) - xopt.getEntry(i));
-                    sl.setEntry(i, sl.getEntry(i) - xopt.getEntry(i));
-                    su.setEntry(i, su.getEntry(i) - xopt.getEntry(i));
-                    xopt.setEntry(i, ZERO);
+                    originShift.setEntry(i, originShift.getEntry(i) + trustRegionCenterOffset.getEntry(i));
+                    newPoint.setEntry(i, newPoint.getEntry(i) - trustRegionCenterOffset.getEntry(i));
+                    lowerDifference.setEntry(i, lowerDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i));
+                    upperDifference.setEntry(i, upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i));
+                    trustRegionCenterOffset.setEntry(i, ZERO);
                 }
                 xoptsq = ZERO;
             }
@@ -652,6 +642,7 @@ public class BOBYQAOptimizer
 
         }
         case 210: {
+            printState(210); // XXX
             // Pick two alternative vectors of variables, relative to XBASE, that
             // are suitable as new positions of the KNEW-th interpolation point.
             // Firstly, XNEW is set to the point on a line through XOPT and another
@@ -663,14 +654,12 @@ public class BOBYQAOptimizer
             // being returned in CAUCHY. The choice between these alternatives is
             // going to be made when the denominator is calculated.
 
-            final double[] alphaCauchy = altmov(xpt, xopt,
-                                                bmat, zmat,
-                                                sl, su, knew, adelt, xnew, xalt);
+            final double[] alphaCauchy = altmov(knew, adelt);
             alpha = alphaCauchy[0];
             cauchy = alphaCauchy[1];
 
             for (int i = 0; i < n; i++) {
-                d__.setEntry(i, xnew.getEntry(i) - xopt.getEntry(i));
+                trialStepPoint.setEntry(i, newPoint.getEntry(i) - trustRegionCenterOffset.getEntry(i));
             }
 
             // Calculate VLAG and BETA for the current choice of D. The scalar
@@ -679,28 +668,29 @@ public class BOBYQAOptimizer
 
         }
         case 230: {
+            printState(230); // XXX
             for (int k = 0; k < npt; k++) {
                 double suma = ZERO;
                 double sumb = ZERO;
                 double sum = ZERO;
                 for (int j = 0; j < n; j++) {
-                    suma += xpt.getEntry(k, j) * d__.getEntry(j);
-                    sumb += xpt.getEntry(k, j) * xopt.getEntry(j);
-                    sum += bmat.getEntry(k, j) * d__.getEntry(j);
+                    suma += interpolationPoints.getEntry(k, j) * trialStepPoint.getEntry(j);
+                    sumb += interpolationPoints.getEntry(k, j) * trustRegionCenterOffset.getEntry(j);
+                    sum += bMatrix.getEntry(k, j) * trialStepPoint.getEntry(j);
                 }
                 work3.setEntry(k, suma * (HALF * suma + sumb));
-                vlag.setEntry(k, sum);
+                lagrangeValuesAtNewPoint.setEntry(k, sum);
                 work2.setEntry(k, suma);
             }
             beta = ZERO;
             for (int m = 0; m < nptm; m++) {
                 double sum = ZERO;
                 for (int k = 0; k < npt; k++) {
-                    sum += zmat.getEntry(k, m) * work3.getEntry(k);
+                    sum += zMatrix.getEntry(k, m) * work3.getEntry(k);
                 }
                 beta -= sum * sum;
                 for (int k = 0; k < npt; k++) {
-                    vlag.setEntry(k, vlag.getEntry(k) + sum * zmat.getEntry(k, m));
+                    lagrangeValuesAtNewPoint.setEntry(k, lagrangeValuesAtNewPoint.getEntry(k) + sum * zMatrix.getEntry(k, m));
                 }
             }
             dsq = ZERO;
@@ -708,28 +698,28 @@ public class BOBYQAOptimizer
             double dx = ZERO;
             for (int j = 0; j < n; j++) {
                 // Computing 2nd power
-                final double d1 = d__.getEntry(j);
+                final double d1 = trialStepPoint.getEntry(j);
                 dsq += d1 * d1;
                 double sum = ZERO;
                 for (int k = 0; k < npt; k++) {
-                    sum += work3.getEntry(k) * bmat.getEntry(k, j);
+                    sum += work3.getEntry(k) * bMatrix.getEntry(k, j);
                 }
-                bsum += sum * d__.getEntry(j);
+                bsum += sum * trialStepPoint.getEntry(j);
                 final int jp = npt + j;
                 for (int i = 0; i < n; i++) {
-                    sum += bmat.getEntry(jp, i) * d__.getEntry(i);
+                    sum += bMatrix.getEntry(jp, i) * trialStepPoint.getEntry(i);
                 }
-                vlag.setEntry(jp, sum);
-                bsum += sum * d__.getEntry(j);
-                dx += d__.getEntry(j) * xopt.getEntry(j);
+                lagrangeValuesAtNewPoint.setEntry(jp, sum);
+                bsum += sum * trialStepPoint.getEntry(j);
+                dx += trialStepPoint.getEntry(j) * trustRegionCenterOffset.getEntry(j);
             }
 
             beta = dx * dx + dsq * (xoptsq + dx + dx + HALF * dsq) + beta - bsum; // Original
             // beta += dx * dx + dsq * (xoptsq + dx + dx + HALF * dsq) - bsum; // XXX "testAckley" and "testDiffPow" fail.
             // beta = dx * dx + dsq * (xoptsq + 2 * dx + HALF * dsq) + beta - bsum; // XXX "testDiffPow" fails.
 
-            vlag.setEntry(trustRegionCenterInterpolationPointIndex,
-                          vlag.getEntry(trustRegionCenterInterpolationPointIndex) + ONE);
+            lagrangeValuesAtNewPoint.setEntry(trustRegionCenterInterpolationPointIndex,
+                          lagrangeValuesAtNewPoint.getEntry(trustRegionCenterInterpolationPointIndex) + ONE);
 
             // If NTRITS is zero, the denominator may be increased by replacing
             // the step D of ALTMOV by a Cauchy step. Then RESCUE may be called if
@@ -737,12 +727,12 @@ public class BOBYQAOptimizer
 
             if (ntrits == 0) {
                 // Computing 2nd power
-                final double d1 = vlag.getEntry(knew);
+                final double d1 = lagrangeValuesAtNewPoint.getEntry(knew);
                 denom = d1 * d1 + alpha * beta;
                 if (denom < cauchy && cauchy > ZERO) {
                     for (int i = 0; i < n; i++) {
-                        xnew.setEntry(i, xalt.getEntry(i));
-                        d__.setEntry(i, xnew.getEntry(i) - xopt.getEntry(i));
+                        newPoint.setEntry(i, alternativeNewPoint.getEntry(i));
+                        trialStepPoint.setEntry(i, newPoint.getEntry(i) - trustRegionCenterOffset.getEntry(i));
                     }
                     cauchy = ZERO; // XXX Useful statement?
                     state = 230; break;
@@ -765,16 +755,16 @@ public class BOBYQAOptimizer
                     double hdiag = ZERO;
                     for (int m = 0; m < nptm; m++) {
                         // Computing 2nd power
-                        final double d1 = zmat.getEntry(k, m);
+                        final double d1 = zMatrix.getEntry(k, m);
                         hdiag += d1 * d1;
                     }
                     // Computing 2nd power
-                    final double d2 = vlag.getEntry(k);
+                    final double d2 = lagrangeValuesAtNewPoint.getEntry(k);
                     final double den = beta * hdiag + d2 * d2;
                     distsq = ZERO;
                     for (int j = 0; j < n; j++) {
                         // Computing 2nd power
-                        final double d3 = xpt.getEntry(k, j) - xopt.getEntry(j);
+                        final double d3 = interpolationPoints.getEntry(k, j) - trustRegionCenterOffset.getEntry(j);
                         distsq += d3 * d3;
                     }
                     // Computing MAX
@@ -788,7 +778,7 @@ public class BOBYQAOptimizer
                     }
                     // Computing MAX
                     // Computing 2nd power
-                    final double d5 = vlag.getEntry(k);
+                    final double d5 = lagrangeValuesAtNewPoint.getEntry(k);
                     biglsq = Math.max(biglsq, temp * (d5 * d5));
                 }
             }
@@ -801,18 +791,19 @@ public class BOBYQAOptimizer
 
         }
         case 360: {
+            printState(360); // XXX
             for (int i = 0; i < n; i++) {
                 // Computing MIN
                 // Computing MAX
                 final double d3 = lowerBound[i];
-                final double d4 = xbase.getEntry(i) + xnew.getEntry(i);
+                final double d4 = originShift.getEntry(i) + newPoint.getEntry(i);
                 final double d1 = Math.max(d3, d4);
                 final double d2 = upperBound[i];
                 currentBest.setEntry(i, Math.min(d1, d2));
-                if (xnew.getEntry(i) == sl.getEntry(i)) {
+                if (newPoint.getEntry(i) == lowerDifference.getEntry(i)) {
                     currentBest.setEntry(i, lowerBound[i]);
                 }
-                if (xnew.getEntry(i) == su.getEntry(i)) {
+                if (newPoint.getEntry(i) == upperDifference.getEntry(i)) {
                     currentBest.setEntry(i, upperBound[i]);
                 }
             }
@@ -829,17 +820,17 @@ public class BOBYQAOptimizer
             // Use the quadratic model to predict the change in F due to the step D,
             //   and set DIFF to the error of this prediction.
 
-            final double fopt = fval.getEntry(trustRegionCenterInterpolationPointIndex);
+            final double fopt = fAtInterpolationPoints.getEntry(trustRegionCenterInterpolationPointIndex);
             double vquad = ZERO;
             int ih = 0;
             for (int j = 0; j < n; j++) {
-                vquad += d__.getEntry(j) * gopt.getEntry(j);
+                vquad += trialStepPoint.getEntry(j) * gradientAtTrustRegionCenter.getEntry(j);
                 for (int i = 0; i <= j; i++) {
-                    double temp = d__.getEntry(i) * d__.getEntry(j);
+                    double temp = trialStepPoint.getEntry(i) * trialStepPoint.getEntry(j);
                     if (i == j) {
                         temp *= HALF;
                     }
-                    vquad += hq.getEntry(ih) * temp;
+                    vquad += modelSecondDerivativesValues.getEntry(ih) * temp;
                     ih++;
                }
             }
@@ -847,7 +838,7 @@ public class BOBYQAOptimizer
                 // Computing 2nd power
                 final double d1 = work2.getEntry(k);
                 final double d2 = d1 * d1; // "d1" must be squared first to prevent test failures.
-                vquad += HALF * pq.getEntry(k) * d2;
+                vquad += HALF * modelSecondDerivativesParameters.getEntry(k) * d2;
             }
             final double diff = f - fopt - vquad;
             diffc = diffb;
@@ -892,16 +883,16 @@ public class BOBYQAOptimizer
                         double hdiag = ZERO;
                         for (int m = 0; m < nptm; m++) {
                             // Computing 2nd power
-                            final double d1 = zmat.getEntry(k, m);
+                            final double d1 = zMatrix.getEntry(k, m);
                             hdiag += d1 * d1;
                         }
                         // Computing 2nd power
-                        final double d1 = vlag.getEntry(k);
+                        final double d1 = lagrangeValuesAtNewPoint.getEntry(k);
                         final double den = beta * hdiag + d1 * d1;
                         distsq = ZERO;
                         for (int j = 0; j < n; j++) {
                             // Computing 2nd power
-                            final double d2 = xpt.getEntry(k, j) - xnew.getEntry(j);
+                            final double d2 = interpolationPoints.getEntry(k, j) - newPoint.getEntry(j);
                             distsq += d2 * d2;
                         }
                         // Computing MAX
@@ -915,7 +906,7 @@ public class BOBYQAOptimizer
                         }
                         // Computing MAX
                         // Computing 2nd power
-                        final double d4 = vlag.getEntry(k);
+                        final double d4 = lagrangeValuesAtNewPoint.getEntry(k);
                         final double d5 = temp * (d4 * d4);
                         biglsq = Math.max(biglsq, d5);
                     }
@@ -929,50 +920,49 @@ public class BOBYQAOptimizer
             // Update BMAT and ZMAT, so that the KNEW-th interpolation point can be
             // moved. Also update the second derivative terms of the model.
 
-            update(bmat, zmat, vlag,
-                   beta, denom, knew);
+            update(beta, denom, knew);
 
             ih = 0;
-            final double pqold = pq.getEntry(knew);
-            pq.setEntry(knew, ZERO);
+            final double pqold = modelSecondDerivativesParameters.getEntry(knew);
+            modelSecondDerivativesParameters.setEntry(knew, ZERO);
             for (int i = 0; i < n; i++) {
-                final double temp = pqold * xpt.getEntry(knew, i);
+                final double temp = pqold * interpolationPoints.getEntry(knew, i);
                 for (int j = 0; j <= i; j++) {
-                    hq.setEntry(ih, hq.getEntry(ih) + temp * xpt.getEntry(knew, j));
+                    modelSecondDerivativesValues.setEntry(ih, modelSecondDerivativesValues.getEntry(ih) + temp * interpolationPoints.getEntry(knew, j));
                     ih++;
                 }
             }
             for (int m = 0; m < nptm; m++) {
-                final double temp = diff * zmat.getEntry(knew, m);
+                final double temp = diff * zMatrix.getEntry(knew, m);
                 for (int k = 0; k < npt; k++) {
-                    pq.setEntry(k, pq.getEntry(k) + temp * zmat.getEntry(k, m));
+                    modelSecondDerivativesParameters.setEntry(k, modelSecondDerivativesParameters.getEntry(k) + temp * zMatrix.getEntry(k, m));
                 }
             }
 
             // Include the new interpolation point, and make the changes to GOPT at
             // the old XOPT that are caused by the updating of the quadratic model.
 
-            fval.setEntry(knew,  f);
+            fAtInterpolationPoints.setEntry(knew,  f);
             for (int i = 0; i < n; i++) {
-                xpt.setEntry(knew, i, xnew.getEntry(i));
-                work1.setEntry(i, bmat.getEntry(knew, i));
+                interpolationPoints.setEntry(knew, i, newPoint.getEntry(i));
+                work1.setEntry(i, bMatrix.getEntry(knew, i));
             }
             for (int k = 0; k < npt; k++) {
                 double suma = ZERO;
                 for (int m = 0; m < nptm; m++) {
-                    suma += zmat.getEntry(knew, m) * zmat.getEntry(k, m);
+                    suma += zMatrix.getEntry(knew, m) * zMatrix.getEntry(k, m);
                 }
                 double sumb = ZERO;
                 for (int j = 0; j < n; j++) {
-                    sumb += xpt.getEntry(k, j) * xopt.getEntry(j);
+                    sumb += interpolationPoints.getEntry(k, j) * trustRegionCenterOffset.getEntry(j);
                 }
                 final double temp = suma * sumb;
                 for (int i = 0; i < n; i++) {
-                    work1.setEntry(i, work1.getEntry(i) + temp * xpt.getEntry(k, i));
+                    work1.setEntry(i, work1.getEntry(i) + temp * interpolationPoints.getEntry(k, i));
                 }
             }
             for (int i = 0; i < n; i++) {
-                gopt.setEntry(i, gopt.getEntry(i) + diff * work1.getEntry(i));
+                gradientAtTrustRegionCenter.setEntry(i, gradientAtTrustRegionCenter.getEntry(i) + diff * work1.getEntry(i));
             }
 
             // Update XOPT, GOPT and KOPT if the new calculated F is less than FOPT.
@@ -982,26 +972,26 @@ public class BOBYQAOptimizer
                 xoptsq = ZERO;
                 ih = 0;
                 for (int j = 0; j < n; j++) {
-                    xopt.setEntry(j, xnew.getEntry(j));
+                    trustRegionCenterOffset.setEntry(j, newPoint.getEntry(j));
                     // Computing 2nd power
-                    final double d1 = xopt.getEntry(j);
+                    final double d1 = trustRegionCenterOffset.getEntry(j);
                     xoptsq += d1 * d1;
                     for (int i = 0; i <= j; i++) {
                         if (i < j) {
-                            gopt.setEntry(j, gopt.getEntry(j) + hq.getEntry(ih) * d__.getEntry(i));
+                            gradientAtTrustRegionCenter.setEntry(j, gradientAtTrustRegionCenter.getEntry(j) + modelSecondDerivativesValues.getEntry(ih) * trialStepPoint.getEntry(i));
                         }
-                        gopt.setEntry(i, gopt.getEntry(i) + hq.getEntry(ih) * d__.getEntry(j));
+                        gradientAtTrustRegionCenter.setEntry(i, gradientAtTrustRegionCenter.getEntry(i) + modelSecondDerivativesValues.getEntry(ih) * trialStepPoint.getEntry(j));
                         ih++;
                     }
                 }
                 for (int k = 0; k < npt; k++) {
                     double temp = ZERO;
                     for (int j = 0; j < n; j++) {
-                        temp += xpt.getEntry(k, j) * d__.getEntry(j);
+                        temp += interpolationPoints.getEntry(k, j) * trialStepPoint.getEntry(j);
                     }
-                    temp *= pq.getEntry(k);
+                    temp *= modelSecondDerivativesParameters.getEntry(k);
                     for (int i = 0; i < n; i++) {
-                        gopt.setEntry(i, gopt.getEntry(i) + temp * xpt.getEntry(k, i));
+                        gradientAtTrustRegionCenter.setEntry(i, gradientAtTrustRegionCenter.getEntry(i) + temp * interpolationPoints.getEntry(k, i));
                     }
                 }
             }
@@ -1012,22 +1002,22 @@ public class BOBYQAOptimizer
 
             if (ntrits > 0) {
                 for (int k = 0; k < npt; k++) {
-                    vlag.setEntry(k, fval.getEntry(k) - fval.getEntry(trustRegionCenterInterpolationPointIndex));
+                    lagrangeValuesAtNewPoint.setEntry(k, fAtInterpolationPoints.getEntry(k) - fAtInterpolationPoints.getEntry(trustRegionCenterInterpolationPointIndex));
                     work3.setEntry(k, ZERO);
                 }
                 for (int j = 0; j < nptm; j++) {
                     double sum = ZERO;
                     for (int k = 0; k < npt; k++) {
-                        sum += zmat.getEntry(k, j) * vlag.getEntry(k);
+                        sum += zMatrix.getEntry(k, j) * lagrangeValuesAtNewPoint.getEntry(k);
                     }
                     for (int k = 0; k < npt; k++) {
-                        work3.setEntry(k, work3.getEntry(k) + sum * zmat.getEntry(k, j));
+                        work3.setEntry(k, work3.getEntry(k) + sum * zMatrix.getEntry(k, j));
                     }
                 }
                 for (int k = 0; k < npt; k++) {
                     double sum = ZERO;
                     for (int j = 0; j < n; j++) {
-                        sum += xpt.getEntry(k, j) * xopt.getEntry(j);
+                        sum += interpolationPoints.getEntry(k, j) * trustRegionCenterOffset.getEntry(j);
                     }
                     work2.setEntry(k, work3.getEntry(k));
                     work3.setEntry(k, sum * work3.getEntry(k));
@@ -1037,32 +1027,32 @@ public class BOBYQAOptimizer
                 for (int i = 0; i < n; i++) {
                     double sum = ZERO;
                     for (int k = 0; k < npt; k++) {
-                        sum += bmat.getEntry(k, i) *
-                            vlag.getEntry(k) + xpt.getEntry(k, i) * work3.getEntry(k);
+                        sum += bMatrix.getEntry(k, i) *
+                            lagrangeValuesAtNewPoint.getEntry(k) + interpolationPoints.getEntry(k, i) * work3.getEntry(k);
                     }
-                    if (xopt.getEntry(i) == sl.getEntry(i)) {
+                    if (trustRegionCenterOffset.getEntry(i) == lowerDifference.getEntry(i)) {
                         // Computing MIN
                         // Computing 2nd power
-                        final double d1 = Math.min(ZERO, gopt.getEntry(i));
+                        final double d1 = Math.min(ZERO, gradientAtTrustRegionCenter.getEntry(i));
                         gqsq += d1 * d1;
                         // Computing 2nd power
                         final double d2 = Math.min(ZERO, sum);
                         gisq += d2 * d2;
-                    } else if (xopt.getEntry(i) == su.getEntry(i)) {
+                    } else if (trustRegionCenterOffset.getEntry(i) == upperDifference.getEntry(i)) {
                         // Computing MAX
                         // Computing 2nd power
-                        final double d1 = Math.max(ZERO, gopt.getEntry(i));
+                        final double d1 = Math.max(ZERO, gradientAtTrustRegionCenter.getEntry(i));
                         gqsq += d1 * d1;
                         // Computing 2nd power
                         final double d2 = Math.max(ZERO, sum);
                         gisq += d2 * d2;
                     } else {
                         // Computing 2nd power
-                        final double d1 = gopt.getEntry(i);
+                        final double d1 = gradientAtTrustRegionCenter.getEntry(i);
                         gqsq += d1 * d1;
                         gisq += sum * sum;
                     }
-                    vlag.setEntry(npt + i, sum);
+                    lagrangeValuesAtNewPoint.setEntry(npt + i, sum);
                 }
 
                 // Test whether to replace the new quadratic model by the least Frobenius
@@ -1075,13 +1065,13 @@ public class BOBYQAOptimizer
                 if (itest >= 3) {
                     for (int i = 0, max = Math.max(npt, nh); i < max; i++) {
                         if (i < n) {
-                            gopt.setEntry(i, vlag.getEntry(npt + i));
+                            gradientAtTrustRegionCenter.setEntry(i, lagrangeValuesAtNewPoint.getEntry(npt + i));
                         }
                         if (i < npt) {
-                            pq.setEntry(i, work2.getEntry(i));
+                            modelSecondDerivativesParameters.setEntry(i, work2.getEntry(i));
                         }
                         if (i < nh) {
-                            hq.setEntry(i, ZERO);
+                            modelSecondDerivativesValues.setEntry(i, ZERO);
                         }
                         itest = 0;
                     }
@@ -1110,12 +1100,13 @@ public class BOBYQAOptimizer
             distsq = Math.max(d1 * d1, d2 * d2);
         }
         case 650: {
+            printState(650); // XXX
             knew = -1;
             for (int k = 0; k < npt; k++) {
                 double sum = ZERO;
                 for (int j = 0; j < n; j++) {
                     // Computing 2nd power
-                    final double d1 = xpt.getEntry(k, j) - xopt.getEntry(j);
+                    final double d1 = interpolationPoints.getEntry(k, j) - trustRegionCenterOffset.getEntry(j);
                     sum += d1 * d1;
                 }
                 if (sum > distsq) {
@@ -1161,6 +1152,7 @@ public class BOBYQAOptimizer
             //   next values of RHO and DELTA.
         }
         case 680: {
+            printState(680); // XXX
             if (rho > stoppingTrustRegionRadius) {
                 delta = HALF * rho;
                 ratio = rho / stoppingTrustRegionRadius;
@@ -1185,23 +1177,24 @@ public class BOBYQAOptimizer
             }
         }
         case 720: {
-            if (fval.getEntry(trustRegionCenterInterpolationPointIndex) <= fsave) {
+            printState(720); // XXX
+            if (fAtInterpolationPoints.getEntry(trustRegionCenterInterpolationPointIndex) <= fsave) {
                 for (int i = 0; i < n; i++) {
                     // Computing MIN
                     // Computing MAX
                     final double d3 = lowerBound[i];
-                    final double d4 = xbase.getEntry(i) + xopt.getEntry(i);
+                    final double d4 = originShift.getEntry(i) + trustRegionCenterOffset.getEntry(i);
                     final double d1 = Math.max(d3, d4);
                     final double d2 = upperBound[i];
                     currentBest.setEntry(i, Math.min(d1, d2));
-                    if (xopt.getEntry(i) == sl.getEntry(i)) {
+                    if (trustRegionCenterOffset.getEntry(i) == lowerDifference.getEntry(i)) {
                         currentBest.setEntry(i, lowerBound[i]);
                     }
-                    if (xopt.getEntry(i) == su.getEntry(i)) {
+                    if (trustRegionCenterOffset.getEntry(i) == upperDifference.getEntry(i)) {
                         currentBest.setEntry(i, upperBound[i]);
                     }
                 }
-                f = fval.getEntry(trustRegionCenterInterpolationPointIndex);
+                f = fAtInterpolationPoints.getEntry(trustRegionCenterInterpolationPointIndex);
             }
             return f;
         }
@@ -1242,30 +1235,14 @@ public class BOBYQAOptimizer
      *
      *     Set the first NPT components of W to the leading elements of the
      *     KNEW-th column of the H matrix.
-     * @param xpt
-     * @param xopt
-     * @param bmat
-     * @param zmat
-     * @param sl
-     * @param su
      * @param knew
      * @param adelt
-     * @param xnew
-     * @param xalt
      */
     private double[] altmov(
-            Array2DRowRealMatrix xpt,
-            ArrayRealVector xopt,
-            Array2DRowRealMatrix bmat,
-            Array2DRowRealMatrix zmat,
-            ArrayRealVector sl,
-            ArrayRealVector su,
             int knew,
-            double adelt,
-            ArrayRealVector xnew,
-            ArrayRealVector xalt
+            double adelt
     ) {
-        // System.out.println("altmov"); // XXX
+        printMethod(); // XXX
 
         final int n = currentBest.getDimension();
         final int npt = numberOfInterpolationPoints;
@@ -1280,9 +1257,9 @@ public class BOBYQAOptimizer
             hcol.setEntry(k, ZERO);
         }
         for (int j = 0, max = npt - n - 1; j < max; j++) {
-            final double tmp = zmat.getEntry(knew, j);
+            final double tmp = zMatrix.getEntry(knew, j);
             for (int k = 0; k < npt; k++) {
-                hcol.setEntry(k, hcol.getEntry(k) + tmp * zmat.getEntry(k, j));
+                hcol.setEntry(k, hcol.getEntry(k) + tmp * zMatrix.getEntry(k, j));
             }
         }
         final double alpha = hcol.getEntry(knew);
@@ -1291,16 +1268,16 @@ public class BOBYQAOptimizer
         // Calculate the gradient of the KNEW-th Lagrange function at XOPT.
 
         for (int i = 0; i < n; i++) {
-            glag.setEntry(i, bmat.getEntry(knew, i));
+            glag.setEntry(i, bMatrix.getEntry(knew, i));
         }
         for (int k = 0; k < npt; k++) {
             double tmp = ZERO;
             for (int j = 0; j < n; j++) {
-                tmp += xpt.getEntry(k, j) * xopt.getEntry(j);
+                tmp += interpolationPoints.getEntry(k, j) * trustRegionCenterOffset.getEntry(j);
             }
             tmp *= hcol.getEntry(k);
             for (int i = 0; i < n; i++) {
-                glag.setEntry(i, glag.getEntry(i) + tmp * xpt.getEntry(k, i));
+                glag.setEntry(i, glag.getEntry(i) + tmp * interpolationPoints.getEntry(k, i));
             }
         }
 
@@ -1322,7 +1299,7 @@ public class BOBYQAOptimizer
             double dderiv = ZERO;
             double distsq = ZERO;
             for (int i = 0; i < n; i++) {
-                final double tmp = xpt.getEntry(k, i) - xopt.getEntry(i);
+                final double tmp = interpolationPoints.getEntry(k, i) - trustRegionCenterOffset.getEntry(i);
                 dderiv += glag.getEntry(i) * tmp;
                 distsq += tmp * tmp;
             }
@@ -1335,27 +1312,27 @@ public class BOBYQAOptimizer
             // Revise SLBD and SUBD if necessary because of the bounds in SL and SU.
 
             for (int i = 0; i < n; i++) {
-                final double tmp = xpt.getEntry(k, i) - xopt.getEntry(i);
+                final double tmp = interpolationPoints.getEntry(k, i) - trustRegionCenterOffset.getEntry(i);
                 if (tmp > ZERO) {
-                    if (slbd * tmp < sl.getEntry(i) - xopt.getEntry(i)) {
-                        slbd = (sl.getEntry(i) - xopt.getEntry(i)) / tmp;
+                    if (slbd * tmp < lowerDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) {
+                        slbd = (lowerDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) / tmp;
                         ilbd = -i - 1;
                     }
-                    if (subd * tmp > su.getEntry(i) - xopt.getEntry(i)) {
+                    if (subd * tmp > upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) {
                         // Computing MAX
                         subd = Math.max(sumin,
-                                        (su.getEntry(i) - xopt.getEntry(i)) / tmp);
+                                        (upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) / tmp);
                         iubd = i + 1;
                     }
                 } else if (tmp < ZERO) {
-                    if (slbd * tmp > su.getEntry(i) - xopt.getEntry(i)) {
-                        slbd = (su.getEntry(i) - xopt.getEntry(i)) / tmp;
+                    if (slbd * tmp > upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) {
+                        slbd = (upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) / tmp;
                         ilbd = i + 1;
                     }
-                    if (subd * tmp < sl.getEntry(i) - xopt.getEntry(i)) {
+                    if (subd * tmp < lowerDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) {
                         // Computing MAX
                         subd = Math.max(sumin,
-                                        (sl.getEntry(i) - xopt.getEntry(i)) / tmp);
+                                        (lowerDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i)) / tmp);
                         iubd = -i - 1;
                     }
                 }
@@ -1423,15 +1400,15 @@ public class BOBYQAOptimizer
         // Construct XNEW in a way that satisfies the bound constraints exactly.
 
         for (int i = 0; i < n; i++) {
-            final double tmp = xopt.getEntry(i) + stpsav * (xpt.getEntry(ksav, i) - xopt.getEntry(i));
-            xnew.setEntry(i, Math.max(sl.getEntry(i),
-                                      Math.min(su.getEntry(i), tmp)));
+            final double tmp = trustRegionCenterOffset.getEntry(i) + stpsav * (interpolationPoints.getEntry(ksav, i) - trustRegionCenterOffset.getEntry(i));
+            newPoint.setEntry(i, Math.max(lowerDifference.getEntry(i),
+                                      Math.min(upperDifference.getEntry(i), tmp)));
         }
         if (ibdsav < 0) {
-            xnew.setEntry(-ibdsav - 1, sl.getEntry(-ibdsav - 1));
+            newPoint.setEntry(-ibdsav - 1, lowerDifference.getEntry(-ibdsav - 1));
         }
         if (ibdsav > 0) {
-            xnew.setEntry(ibdsav - 1, su.getEntry(ibdsav - 1));
+            newPoint.setEntry(ibdsav - 1, upperDifference.getEntry(ibdsav - 1));
         }
 
         // Prepare for the iterative method that assembles the constrained Cauchy
@@ -1448,8 +1425,8 @@ public class BOBYQAOptimizer
             for (int i = 0; i < n; i++) {
                 final double glagValue = glag.getEntry(i);
                 work1.setEntry(i, ZERO);
-                if (Math.min(xopt.getEntry(i) - sl.getEntry(i), glagValue) > ZERO ||
-                    Math.max(xopt.getEntry(i) - su.getEntry(i), glagValue) < ZERO) {
+                if (Math.min(trustRegionCenterOffset.getEntry(i) - lowerDifference.getEntry(i), glagValue) > ZERO ||
+                    Math.max(trustRegionCenterOffset.getEntry(i) - upperDifference.getEntry(i), glagValue) < ZERO) {
                     work1.setEntry(i, bigstp);
                     // Computing 2nd power
                     ggfree += glagValue * glagValue;
@@ -1467,14 +1444,14 @@ public class BOBYQAOptimizer
                 ggfree = ZERO;
                 for (int i = 0; i < n; i++) {
                     if (work1.getEntry(i) == bigstp) {
-                        final double tmp2 = xopt.getEntry(i) - step * glag.getEntry(i);
-                        if (tmp2 <= sl.getEntry(i)) {
-                            work1.setEntry(i, sl.getEntry(i) - xopt.getEntry(i));
+                        final double tmp2 = trustRegionCenterOffset.getEntry(i) - step * glag.getEntry(i);
+                        if (tmp2 <= lowerDifference.getEntry(i)) {
+                            work1.setEntry(i, lowerDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i));
                             // Computing 2nd power
                             final double d1 = work1.getEntry(i);
                             wfixsq += d1 * d1;
-                        } else if (tmp2 >= su.getEntry(i)) {
-                            work1.setEntry(i, su.getEntry(i) - xopt.getEntry(i));
+                        } else if (tmp2 >= upperDifference.getEntry(i)) {
+                            work1.setEntry(i, upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i));
                             // Computing 2nd power
                             final double d1 = work1.getEntry(i);
                             wfixsq += d1 * d1;
@@ -1495,15 +1472,15 @@ public class BOBYQAOptimizer
                 final double glagValue = glag.getEntry(i);
                 if (work1.getEntry(i) == bigstp) {
                     work1.setEntry(i, -step * glagValue);
-                    final double min = Math.min(su.getEntry(i),
-                                                xopt.getEntry(i) + work1.getEntry(i));
-                    xalt.setEntry(i, Math.max(sl.getEntry(i), min));
+                    final double min = Math.min(upperDifference.getEntry(i),
+                                                trustRegionCenterOffset.getEntry(i) + work1.getEntry(i));
+                    alternativeNewPoint.setEntry(i, Math.max(lowerDifference.getEntry(i), min));
                 } else if (work1.getEntry(i) == ZERO) {
-                    xalt.setEntry(i, xopt.getEntry(i));
+                    alternativeNewPoint.setEntry(i, trustRegionCenterOffset.getEntry(i));
                 } else if (glagValue > ZERO) {
-                    xalt.setEntry(i, sl.getEntry(i));
+                    alternativeNewPoint.setEntry(i, lowerDifference.getEntry(i));
                 } else {
-                    xalt.setEntry(i, su.getEntry(i));
+                    alternativeNewPoint.setEntry(i, upperDifference.getEntry(i));
                 }
                 gw += glagValue * work1.getEntry(i);
             }
@@ -1517,7 +1494,7 @@ public class BOBYQAOptimizer
             for (int k = 0; k < npt; k++) {
                 double tmp = ZERO;
                 for (int j = 0; j < n; j++) {
-                    tmp += xpt.getEntry(k, j) * work1.getEntry(j);
+                    tmp += interpolationPoints.getEntry(k, j) * work1.getEntry(j);
                 }
                 curv += hcol.getEntry(k) * tmp * tmp;
             }
@@ -1528,9 +1505,9 @@ public class BOBYQAOptimizer
                 curv < -gw * (ONE + Math.sqrt(TWO))) {
                 final double scale = -gw / curv;
                 for (int i = 0; i < n; i++) {
-                    final double tmp = xopt.getEntry(i) + scale * work1.getEntry(i);
-                    xalt.setEntry(i, Math.max(sl.getEntry(i),
-                                              Math.min(su.getEntry(i), tmp)));
+                    final double tmp = trustRegionCenterOffset.getEntry(i) + scale * work1.getEntry(i);
+                    alternativeNewPoint.setEntry(i, Math.max(lowerDifference.getEntry(i),
+                                              Math.min(upperDifference.getEntry(i), tmp)));
                 }
                 // Computing 2nd power
                 final double d1 = HALF * gw * scale;
@@ -1548,7 +1525,7 @@ public class BOBYQAOptimizer
             if (iflag == 0) {
                 for (int i = 0; i < n; i++) {
                     glag.setEntry(i, -glag.getEntry(i));
-                    work2.setEntry(i, xalt.getEntry(i));
+                    work2.setEntry(i, alternativeNewPoint.getEntry(i));
                 }
                 csave = cauchy;
                 iflag = 1;
@@ -1558,7 +1535,7 @@ public class BOBYQAOptimizer
         }
         if (csave > cauchy) {
             for (int i = 0; i < n; i++) {
-                xalt.setEntry(i, work2.getEntry(i));
+                alternativeNewPoint.setEntry(i, work2.getEntry(i));
             }
             cauchy = csave;
         }
@@ -1585,36 +1562,13 @@ public class BOBYQAOptimizer
      *     KOPT will be such that the least calculated value of F so far is at
      *       the point XPT(KOPT,.)+XBASE in the space of the variables.
      *
-     * @param currentBest
-     * @param xbase
-     * @param xpt
-     * @param fval
-     * @param gopt
-     * @param hq
-     * @param pq
-     * @param bmat
-     * @param zmat
-     * @param sl
-     * @param su
-     */
-    private void prelim(
-            ArrayRealVector currentBest,
-            ArrayRealVector xbase,
-            Array2DRowRealMatrix xpt,
-            ArrayRealVector fval,
-            ArrayRealVector gopt,
-            ArrayRealVector hq,
-            ArrayRealVector pq,
-            Array2DRowRealMatrix bmat,
-            Array2DRowRealMatrix zmat,
-            ArrayRealVector sl,
-            ArrayRealVector su
-    ) {
-        // System.out.println("prelim"); // XXX
+     */
+    private void prelim() {
+        printMethod(); // XXX
 
         final int n = currentBest.getDimension();
         final int npt = numberOfInterpolationPoints;
-        final int ndim = bmat.getRowDimension();
+        final int ndim = bMatrix.getRowDimension();
 
         final double rhosq = initialTrustRegionRadius * initialTrustRegionRadius;
         final double recip = 1d / rhosq;
@@ -1624,21 +1578,21 @@ public class BOBYQAOptimizer
         // elements of XPT, BMAT, HQ, PQ and ZMAT to zero.
 
         for (int j = 0; j < n; j++) {
-            xbase.setEntry(j, currentBest.getEntry(j));
+            originShift.setEntry(j, currentBest.getEntry(j));
             for (int k = 0; k < npt; k++) {
-                xpt.setEntry(k, j, ZERO);
+                interpolationPoints.setEntry(k, j, ZERO);
             }
             for (int i = 0; i < ndim; i++) {
-                bmat.setEntry(i, j, ZERO);
+                bMatrix.setEntry(i, j, ZERO);
             }
         }
         for (int i = 0, max = n * np / 2; i < max; i++) {
-            hq.setEntry(i, ZERO);
+            modelSecondDerivativesValues.setEntry(i, ZERO);
         }
         for (int k = 0; k < npt; k++) {
-            pq.setEntry(k, ZERO);
+            modelSecondDerivativesParameters.setEntry(k, ZERO);
             for (int j = 0, max = npt - np; j < max; j++) {
-                zmat.setEntry(k, j, ZERO);
+                zMatrix.setEntry(k, j, ZERO);
             }
         }
 
@@ -1660,23 +1614,23 @@ public class BOBYQAOptimizer
                 if (nfm >= 1 &&
                     nfm <= n) {
                     stepa = initialTrustRegionRadius;
-                    if (su.getEntry(nfmm) == ZERO) {
+                    if (upperDifference.getEntry(nfmm) == ZERO) {
                         stepa = -stepa;
                         throw new PathIsExploredException(); // XXX
                     }
-                    xpt.setEntry(nfm, nfmm, stepa);
+                    interpolationPoints.setEntry(nfm, nfmm, stepa);
                 } else if (nfm > n) {
-                    stepa = xpt.getEntry(nfx, nfxm);
+                    stepa = interpolationPoints.getEntry(nfx, nfxm);
                     stepb = -initialTrustRegionRadius;
-                    if (sl.getEntry(nfxm) == ZERO) {
-                        stepb = Math.min(TWO * initialTrustRegionRadius, su.getEntry(nfxm));
+                    if (lowerDifference.getEntry(nfxm) == ZERO) {
+                        stepb = Math.min(TWO * initialTrustRegionRadius, upperDifference.getEntry(nfxm));
                         throw new PathIsExploredException(); // XXX
                     }
-                    if (su.getEntry(nfxm) == ZERO) {
-                        stepb = Math.max(-TWO * initialTrustRegionRadius, sl.getEntry(nfxm));
+                    if (upperDifference.getEntry(nfxm) == ZERO) {
+                        stepb = Math.max(-TWO * initialTrustRegionRadius, lowerDifference.getEntry(nfxm));
                         throw new PathIsExploredException(); // XXX
                     }
-                    xpt.setEntry(nfm, nfxm, stepb);
+                    interpolationPoints.setEntry(nfm, nfxm, stepb);
                 }
             } else {
                 final int tmp1 = (nfm - np) / n;
@@ -1688,8 +1642,8 @@ public class BOBYQAOptimizer
                     ipt = tmp2;
                     throw new PathIsExploredException(); // XXX
                 }
-                xpt.setEntry(nfm, ipt, xpt.getEntry(ipt, ipt));
-                xpt.setEntry(nfm, jpt, xpt.getEntry(jpt, jpt));
+                interpolationPoints.setEntry(nfm, ipt, interpolationPoints.getEntry(ipt, ipt));
+                interpolationPoints.setEntry(nfm, jpt, interpolationPoints.getEntry(jpt, jpt));
             }
 
             // Calculate the next value of F. The least function value so far and
@@ -1697,12 +1651,12 @@ public class BOBYQAOptimizer
 
             for (int j = 0; j < n; j++) {
                 currentBest.setEntry(j, Math.min(Math.max(lowerBound[j],
-                                                          xbase.getEntry(j) + xpt.getEntry(nfm, j)),
+                                                          originShift.getEntry(j) + interpolationPoints.getEntry(nfm, j)),
                                                  upperBound[j]));
-                if (xpt.getEntry(nfm, j) == sl.getEntry(j)) {
+                if (interpolationPoints.getEntry(nfm, j) == lowerDifference.getEntry(j)) {
                     currentBest.setEntry(j, lowerBound[j]);
                 }
-                if (xpt.getEntry(nfm, j) == su.getEntry(j)) {
+                if (interpolationPoints.getEntry(nfm, j) == upperDifference.getEntry(j)) {
                     currentBest.setEntry(j, upperBound[j]);
                 }
             }
@@ -1710,12 +1664,12 @@ public class BOBYQAOptimizer
             final double objectiveValue = computeObjectiveValue(currentBest.toArray());
             final double f = isMinimize ? objectiveValue : -objectiveValue;
             final int numEval = getEvaluations(); // nfm + 1
-            fval.setEntry(nfm, f);
+            fAtInterpolationPoints.setEntry(nfm, f);
 
             if (numEval == 1) {
                 fbeg = f;
                 trustRegionCenterInterpolationPointIndex = 0;
-            } else if (f < fval.getEntry(trustRegionCenterInterpolationPointIndex)) {
+            } else if (f < fAtInterpolationPoints.getEntry(trustRegionCenterInterpolationPointIndex)) {
                 trustRegionCenterInterpolationPointIndex = nfm;
             }
 
@@ -1728,54 +1682,54 @@ public class BOBYQAOptimizer
             if (numEval <= 2 * n + 1) {
                 if (numEval >= 2 &&
                     numEval <= n + 1) {
-                    gopt.setEntry(nfmm, (f - fbeg) / stepa);
+                    gradientAtTrustRegionCenter.setEntry(nfmm, (f - fbeg) / stepa);
                     if (npt < numEval + n) {
                         final double oneOverStepA = ONE / stepa;
-                        bmat.setEntry(0, nfmm, -oneOverStepA);
-                        bmat.setEntry(nfm, nfmm, oneOverStepA);
-                        bmat.setEntry(npt + nfmm, nfmm, -HALF * rhosq);
+                        bMatrix.setEntry(0, nfmm, -oneOverStepA);
+                        bMatrix.setEntry(nfm, nfmm, oneOverStepA);
+                        bMatrix.setEntry(npt + nfmm, nfmm, -HALF * rhosq);
                         throw new PathIsExploredException(); // XXX
                     }
                 } else if (numEval >= n + 2) {
                     final int ih = nfx * (nfx + 1) / 2 - 1;
                     final double tmp = (f - fbeg) / stepb;
                     final double diff = stepb - stepa;
-                    hq.setEntry(ih, TWO * (tmp - gopt.getEntry(nfxm)) / diff);
-                    gopt.setEntry(nfxm, (gopt.getEntry(nfxm) * stepb - tmp * stepa) / diff);
+                    modelSecondDerivativesValues.setEntry(ih, TWO * (tmp - gradientAtTrustRegionCenter.getEntry(nfxm)) / diff);
+                    gradientAtTrustRegionCenter.setEntry(nfxm, (gradientAtTrustRegionCenter.getEntry(nfxm) * stepb - tmp * stepa) / diff);
                     if (stepa * stepb < ZERO) {
-                        if (f < fval.getEntry(nfm - n)) {
-                            fval.setEntry(nfm, fval.getEntry(nfm - n));
-                            fval.setEntry(nfm - n, f);
+                        if (f < fAtInterpolationPoints.getEntry(nfm - n)) {
+                            fAtInterpolationPoints.setEntry(nfm, fAtInterpolationPoints.getEntry(nfm - n));
+                            fAtInterpolationPoints.setEntry(nfm - n, f);
                             if (trustRegionCenterInterpolationPointIndex == nfm) {
                                 trustRegionCenterInterpolationPointIndex = nfm - n;
                             }
-                            xpt.setEntry(nfm - n, nfxm, stepb);
-                            xpt.setEntry(nfm, nfxm, stepa);
+                            interpolationPoints.setEntry(nfm - n, nfxm, stepb);
+                            interpolationPoints.setEntry(nfm, nfxm, stepa);
                         }
                     }
-                    bmat.setEntry(0, nfxm, -(stepa + stepb) / (stepa * stepb));
-                    bmat.setEntry(nfm, nfxm, -HALF / xpt.getEntry(nfm - n, nfxm));
-                    bmat.setEntry(nfm - n, nfxm,
-                                  -bmat.getEntry(0, nfxm) - bmat.getEntry(nfm, nfxm));
-                    zmat.setEntry(0, nfxm, Math.sqrt(TWO) / (stepa * stepb));
-                    zmat.setEntry(nfm, nfxm, Math.sqrt(HALF) / rhosq);
-                    // zmat.setEntry(nfm, nfxm, Math.sqrt(HALF) * recip); // XXX "testAckley" and "testDiffPow" fail.
-                    zmat.setEntry(nfm - n, nfxm,
-                                  -zmat.getEntry(0, nfxm) - zmat.getEntry(nfm, nfxm));
+                    bMatrix.setEntry(0, nfxm, -(stepa + stepb) / (stepa * stepb));
+                    bMatrix.setEntry(nfm, nfxm, -HALF / interpolationPoints.getEntry(nfm - n, nfxm));
+                    bMatrix.setEntry(nfm - n, nfxm,
+                                  -bMatrix.getEntry(0, nfxm) - bMatrix.getEntry(nfm, nfxm));
+                    zMatrix.setEntry(0, nfxm, Math.sqrt(TWO) / (stepa * stepb));
+                    zMatrix.setEntry(nfm, nfxm, Math.sqrt(HALF) / rhosq);
+                    // zMatrix.setEntry(nfm, nfxm, Math.sqrt(HALF) * recip); // XXX "testAckley" and "testDiffPow" fail.
+                    zMatrix.setEntry(nfm - n, nfxm,
+                                  -zMatrix.getEntry(0, nfxm) - zMatrix.getEntry(nfm, nfxm));
                 }
 
                 // Set the off-diagonal second derivatives of the Lagrange functions and
                 // the initial quadratic model.
 
             } else {
-                zmat.setEntry(0, nfxm, recip);
-                zmat.setEntry(nfm, nfxm, recip);
-                zmat.setEntry(ipt, nfxm, -recip);
-                zmat.setEntry(jpt, nfxm, -recip);
+                zMatrix.setEntry(0, nfxm, recip);
+                zMatrix.setEntry(nfm, nfxm, recip);
+                zMatrix.setEntry(ipt, nfxm, -recip);
+                zMatrix.setEntry(jpt, nfxm, -recip);
 
                 final int ih = ipt * (ipt - 1) / 2 + jpt - 1;
-                final double tmp = xpt.getEntry(nfm, ipt - 1) * xpt.getEntry(nfm, jpt - 1);
-                hq.setEntry(ih, (fbeg - fval.getEntry(ipt) - fval.getEntry(jpt) + f) / tmp);
+                final double tmp = interpolationPoints.getEntry(nfm, ipt - 1) * interpolationPoints.getEntry(nfm, jpt - 1);
+                modelSecondDerivativesValues.setEntry(ih, (fbeg - fAtInterpolationPoints.getEntry(ipt) - fAtInterpolationPoints.getEntry(jpt) + f) / tmp);
                 throw new PathIsExploredException(); // XXX
             }
         } while (getEvaluations() < npt);
@@ -1821,16 +1775,7 @@ public class BOBYQAOptimizer
      *       gradient searches that are not restricted by any constraints. The
      *       value CRVMIN=-1.0D0 is set, however, if all of these searches are
      *       constrained.
-     * @param xpt
-     * @param xopt
-     * @param gopt
-     * @param hq
-     * @param pq
-     * @param sl
-     * @param su
      * @param delta
-     * @param xnew
-     * @param d__
      * @param gnew
      * @param xbdi
      * @param s
@@ -1838,23 +1783,14 @@ public class BOBYQAOptimizer
      * @param hred
      */
     private double[] trsbox(
-            Array2DRowRealMatrix xpt,
-            ArrayRealVector xopt,
-            ArrayRealVector gopt,
-            ArrayRealVector hq,
-            ArrayRealVector pq,
-            ArrayRealVector sl,
-            ArrayRealVector su,
             double delta,
-            ArrayRealVector xnew,
-            ArrayRealVector d__,
             ArrayRealVector gnew,
             ArrayRealVector xbdi,
             ArrayRealVector s,
             ArrayRealVector hs,
             ArrayRealVector hred
     ) {
-        // System.out.println("trsbox"); // XXX
+        printMethod(); // XXX
 
         final int n = currentBest.getDimension();
         final int npt = numberOfInterpolationPoints;
@@ -1893,20 +1829,20 @@ public class BOBYQAOptimizer
         nact = 0;
         for (int i = 0; i < n; i++) {
             xbdi.setEntry(i, ZERO);
-            if (xopt.getEntry(i) <= sl.getEntry(i)) {
-                if (gopt.getEntry(i) >= ZERO) {
+            if (trustRegionCenterOffset.getEntry(i) <= lowerDifference.getEntry(i)) {
+                if (gradientAtTrustRegionCenter.getEntry(i) >= ZERO) {
                     xbdi.setEntry(i, MINUS_ONE);
                 }
-            } else if (xopt.getEntry(i) >= su.getEntry(i)) {
-                if (gopt.getEntry(i) <= ZERO) {
+            } else if (trustRegionCenterOffset.getEntry(i) >= upperDifference.getEntry(i)) {
+                if (gradientAtTrustRegionCenter.getEntry(i) <= ZERO) {
                     xbdi.setEntry(i, ONE);
                 }
             }
             if (xbdi.getEntry(i) != ZERO) {
                 ++nact;
             }
-            d__.setEntry(i, ZERO);
-            gnew.setEntry(i, gopt.getEntry(i));
+            trialStepPoint.setEntry(i, ZERO);
+            gnew.setEntry(i, gradientAtTrustRegionCenter.getEntry(i));
         }
         delsq = delta * delta;
         qred = ZERO;
@@ -1920,12 +1856,13 @@ public class BOBYQAOptimizer
 
         int state = 20;
         for(;;) {
-            // System.out.println("loop in trsbox: state=" + state); // XXX
             switch (state) {
         case 20: {
+            printState(20); // XXX
             beta = ZERO;
         }
         case 30: {
+            printState(30); // XXX
             stepsq = ZERO;
             for (int i = 0; i < n; i++) {
                 if (xbdi.getEntry(i) != ZERO) {
@@ -1958,15 +1895,16 @@ public class BOBYQAOptimizer
             state = 210; break;
         }
         case 50: {
+            printState(50); // XXX
             resid = delsq;
             ds = ZERO;
             shs = ZERO;
             for (int i = 0; i < n; i++) {
                 if (xbdi.getEntry(i) == ZERO) {
                     // Computing 2nd power
-                    final double d1 = d__.getEntry(i);
+                    final double d1 = trialStepPoint.getEntry(i);
                     resid -= d1 * d1;
-                    ds += s.getEntry(i) * d__.getEntry(i);
+                    ds += s.getEntry(i) * trialStepPoint.getEntry(i);
                     shs += s.getEntry(i) * hs.getEntry(i);
                 }
             }
@@ -1991,11 +1929,11 @@ public class BOBYQAOptimizer
             iact = -1;
             for (int i = 0; i < n; i++) {
                 if (s.getEntry(i) != ZERO) {
-                    xsum = xopt.getEntry(i) + d__.getEntry(i);
+                    xsum = trustRegionCenterOffset.getEntry(i) + trialStepPoint.getEntry(i);
                     if (s.getEntry(i) > ZERO) {
-                        temp = (su.getEntry(i) - xsum) / s.getEntry(i);
+                        temp = (upperDifference.getEntry(i) - xsum) / s.getEntry(i);
                     } else {
-                        temp = (sl.getEntry(i) - xsum) / s.getEntry(i);
+                        temp = (lowerDifference.getEntry(i) - xsum) / s.getEntry(i);
                     }
                     if (temp < stplen) {
                         stplen = temp;
@@ -2025,7 +1963,7 @@ public class BOBYQAOptimizer
                         final double d1 = gnew.getEntry(i);
                         gredsq += d1 * d1;
                     }
-                    d__.setEntry(i, d__.getEntry(i) + stplen * s.getEntry(i));
+                    trialStepPoint.setEntry(i, trialStepPoint.getEntry(i) + stplen * s.getEntry(i));
                 }
                 // Computing MAX
                 final double d1 = stplen * (ggsav - HALF * stplen * shs);
@@ -2042,7 +1980,7 @@ public class BOBYQAOptimizer
                     xbdi.setEntry(iact, MINUS_ONE);
                 }
                 // Computing 2nd power
-                final double d1 = d__.getEntry(iact);
+                final double d1 = trialStepPoint.getEntry(iact);
                 delsq -= d1 * d1;
                 if (delsq <= ZERO) {
                     state = 190; break;
@@ -2065,6 +2003,7 @@ public class BOBYQAOptimizer
             }
         }
         case 90: {
+            printState(90); // XXX
             crvmin = ZERO;
 
             // Prepare for the alternative iteration by calculating some scalars
@@ -2073,6 +2012,7 @@ public class BOBYQAOptimizer
 
         }
         case 100: {
+            printState(100); // XXX
             if (nact >= n - 1) {
                 state = 190; break;
             }
@@ -2082,13 +2022,13 @@ public class BOBYQAOptimizer
             for (int i = 0; i < n; i++) {
                 if (xbdi.getEntry(i) == ZERO) {
                     // Computing 2nd power
-                    double d1 = d__.getEntry(i);
+                    double d1 = trialStepPoint.getEntry(i);
                     dredsq += d1 * d1;
-                    dredg += d__.getEntry(i) * gnew.getEntry(i);
+                    dredg += trialStepPoint.getEntry(i) * gnew.getEntry(i);
                     // Computing 2nd power
                     d1 = gnew.getEntry(i);
                     gredsq += d1 * d1;
-                    s.setEntry(i, d__.getEntry(i));
+                    s.setEntry(i, trialStepPoint.getEntry(i));
                 } else {
                     s.setEntry(i, ZERO);
                 }
@@ -2099,6 +2039,7 @@ public class BOBYQAOptimizer
             // and the reduced G that is orthogonal to the reduced D.
         }
         case 120: {
+            printState(120); // XXX
             ++iterc;
             temp = gredsq * dredsq - dredg * dredg;
             if (temp <= qred * 1e-4 * qred) {
@@ -2107,7 +2048,7 @@ public class BOBYQAOptimizer
             temp = Math.sqrt(temp);
             for (int i = 0; i < n; i++) {
                 if (xbdi.getEntry(i) == ZERO) {
-                    s.setEntry(i, (dredg * d__.getEntry(i) - dredsq * gnew.getEntry(i)) / temp);
+                    s.setEntry(i, (dredg * trialStepPoint.getEntry(i) - dredsq * gnew.getEntry(i)) / temp);
                 } else {
                     s.setEntry(i, ZERO);
                 }
@@ -2123,8 +2064,8 @@ public class BOBYQAOptimizer
             iact = -1;
             for (int i = 0; i < n; i++) {
                 if (xbdi.getEntry(i) == ZERO) {
-                    tempa = xopt.getEntry(i) + d__.getEntry(i) - sl.getEntry(i);
-                    tempb = su.getEntry(i) - xopt.getEntry(i) - d__.getEntry(i);
+                    tempa = trustRegionCenterOffset.getEntry(i) + trialStepPoint.getEntry(i) - lowerDifference.getEntry(i);
+                    tempb = upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i) - trialStepPoint.getEntry(i);
                     if (tempa <= ZERO) {
                         ++nact;
                         xbdi.setEntry(i, MINUS_ONE);
@@ -2135,12 +2076,12 @@ public class BOBYQAOptimizer
                         state = 100; break;
                     }
                     // Computing 2nd power
-                    double d1 = d__.getEntry(i);
+                    double d1 = trialStepPoint.getEntry(i);
                     // Computing 2nd power
                     double d2 = s.getEntry(i);
                     ssq = d1 * d1 + d2 * d2;
                     // Computing 2nd power
-                    d1 = xopt.getEntry(i) - sl.getEntry(i);
+                    d1 = trustRegionCenterOffset.getEntry(i) - lowerDifference.getEntry(i);
                     temp = ssq - d1 * d1;
                     if (temp > ZERO) {
                         temp = Math.sqrt(temp) - s.getEntry(i);
@@ -2151,7 +2092,7 @@ public class BOBYQAOptimizer
                         }
                     }
                     // Computing 2nd power
-                    d1 = su.getEntry(i) - xopt.getEntry(i);
+                    d1 = upperDifference.getEntry(i) - trustRegionCenterOffset.getEntry(i);
                     temp = ssq - d1 * d1;
                     if (temp > ZERO) {
                         temp = Math.sqrt(temp) + s.getEntry(i);
@@ -2169,14 +2110,15 @@ public class BOBYQAOptimizer
             state = 210; break;
         }
         case 150: {
+            printState(150); // XXX
             shs = ZERO;
             dhs = ZERO;
             dhd = ZERO;
             for (int i = 0; i < n; i++) {
                 if (xbdi.getEntry(i) == ZERO) {
                     shs += s.getEntry(i) * hs.getEntry(i);
-                    dhs += d__.getEntry(i) * hs.getEntry(i);
-                    dhd += d__.getEntry(i) * hred.getEntry(i);
+                    dhs += trialStepPoint.getEntry(i) * hs.getEntry(i);
+                    dhd += trialStepPoint.getEntry(i) * hred.getEntry(i);
                 }
             }
 
@@ -2230,8 +2172,8 @@ public class BOBYQAOptimizer
             for (int i = 0; i < n; i++) {
                 gnew.setEntry(i, gnew.getEntry(i) + (cth - ONE) * hred.getEntry(i) + sth * hs.getEntry(i));
                 if (xbdi.getEntry(i) == ZERO) {
-                    d__.setEntry(i, cth * d__.getEntry(i) + sth * s.getEntry(i));
-                    dredg += d__.getEntry(i) * gnew.getEntry(i);
+                    trialStepPoint.setEntry(i, cth * trialStepPoint.getEntry(i) + sth * s.getEntry(i));
+                    dredg += trialStepPoint.getEntry(i) * gnew.getEntry(i);
                     // Computing 2nd power
                     final double d1 = gnew.getEntry(i);
                     gredsq += d1 * d1;
@@ -2253,22 +2195,23 @@ public class BOBYQAOptimizer
             }
         }
         case 190: {
+            printState(190); // XXX
             dsq = ZERO;
             for (int i = 0; i < n; i++) {
                 // Computing MAX
                 // Computing MIN
-                final double min = Math.min(xopt.getEntry(i) + d__.getEntry(i),
-                                            su.getEntry(i));
-                xnew.setEntry(i, Math.max(min, sl.getEntry(i)));
+                final double min = Math.min(trustRegionCenterOffset.getEntry(i) + trialStepPoint.getEntry(i),
+                                            upperDifference.getEntry(i));
+                newPoint.setEntry(i, Math.max(min, lowerDifference.getEntry(i)));
                 if (xbdi.getEntry(i) == MINUS_ONE) {
-                    xnew.setEntry(i, sl.getEntry(i));
+                    newPoint.setEntry(i, lowerDifference.getEntry(i));
                 }
                 if (xbdi.getEntry(i) == ONE) {
-                    xnew.setEntry(i, su.getEntry(i));
+                    newPoint.setEntry(i, upperDifference.getEntry(i));
                 }
-                d__.setEntry(i, xnew.getEntry(i) - xopt.getEntry(i));
+                trialStepPoint.setEntry(i, newPoint.getEntry(i) - trustRegionCenterOffset.getEntry(i));
                 // Computing 2nd power
-                final double d1 = d__.getEntry(i);
+                final double d1 = trialStepPoint.getEntry(i);
                 dsq += d1 * d1;
             }
             return new double[] { dsq, crvmin };
@@ -2278,22 +2221,23 @@ public class BOBYQAOptimizer
             // they can be regarded as an external subroutine.
         }
         case 210: {
+            printState(210); // XXX
             int ih = 0;
             for (int j = 0; j < n; j++) {
                 hs.setEntry(j, ZERO);
                 for (int i = 0; i <= j; i++) {
                     if (i < j) {
-                        hs.setEntry(j, hs.getEntry(j) + hq.getEntry(ih) * s.getEntry(i));
+                        hs.setEntry(j, hs.getEntry(j) + modelSecondDerivativesValues.getEntry(ih) * s.getEntry(i));
                     }
-                    hs.setEntry(i, hs.getEntry(i) + hq.getEntry(ih) * s.getEntry(j));
+                    hs.setEntry(i, hs.getEntry(i) + modelSecondDerivativesValues.getEntry(ih) * s.getEntry(j));
                     ih++;
                 }
             }
-            final RealVector tmp = xpt.operate(s).ebeMultiply(pq);
+            final RealVector tmp = interpolationPoints.operate(s).ebeMultiply(modelSecondDerivativesParameters);
             for (int k = 0; k < npt; k++) {
-                if (pq.getEntry(k) != ZERO) {
+                if (modelSecondDerivativesParameters.getEntry(k) != ZERO) {
                     for (int i = 0; i < n; i++) {
-                        hs.setEntry(i, hs.getEntry(i) + tmp.getEntry(k) * xpt.getEntry(k, i));
+                        hs.setEntry(i, hs.getEntry(i) + tmp.getEntry(k) * interpolationPoints.getEntry(k, i));
                     }
                 }
             }
@@ -2325,22 +2269,16 @@ public class BOBYQAOptimizer
      *     with that name, and DENOM is set to the denominator of the updating
      *     formula. Elements of ZMAT may be treated as zero if their moduli are
      *     at most ZTEST. The first NDIM elements of W are used for working space.
-     * @param bmat
-     * @param zmat
-     * @param vlag
      * @param beta
      * @param denom
      * @param knew
      */
     private void update(
-            Array2DRowRealMatrix bmat,
-            Array2DRowRealMatrix zmat,
-            ArrayRealVector vlag,
             double beta,
             double denom,
             int knew
     ) {
-        // System.out.println("update"); // XXX
+        printMethod(); // XXX
 
         final int n = currentBest.getDimension();
         final int npt = numberOfInterpolationPoints;
@@ -2353,7 +2291,7 @@ public class BOBYQAOptimizer
         for (int k = 0; k < npt; k++) {
             for (int j = 0; j < nptm; j++) {
                 // Computing MAX
-                ztest = Math.max(ztest, Math.abs(zmat.getEntry(k, j)));
+                ztest = Math.max(ztest, Math.abs(zMatrix.getEntry(k, j)));
             }
         }
         ztest *= 1e-20;
@@ -2361,56 +2299,56 @@ public class BOBYQAOptimizer
         // Apply the rotations that put zeros in the KNEW-th row of ZMAT.
 
         for (int j = 1; j < nptm; j++) {
-            final double d1 = zmat.getEntry(knew, j);
+            final double d1 = zMatrix.getEntry(knew, j);
             if (Math.abs(d1) > ztest) {
                 // Computing 2nd power
-                final double d2 = zmat.getEntry(knew, 0);
+                final double d2 = zMatrix.getEntry(knew, 0);
                 // Computing 2nd power
-                final double d3 = zmat.getEntry(knew, j);
+                final double d3 = zMatrix.getEntry(knew, j);
                 final double d4 = Math.sqrt(d2 * d2 + d3 * d3);
-                final double d5 = zmat.getEntry(knew, 0) / d4;
-                final double d6 = zmat.getEntry(knew, j) / d4;
+                final double d5 = zMatrix.getEntry(knew, 0) / d4;
+                final double d6 = zMatrix.getEntry(knew, j) / d4;
                 for (int i = 0; i < npt; i++) {
-                    final double d7 = d5 * zmat.getEntry(i, 0) + d6 * zmat.getEntry(i, j);
-                    zmat.setEntry(i, j, d5 * zmat.getEntry(i, j) - d6 * zmat.getEntry(i, 0));
-                    zmat.setEntry(i, 0, d7);
+                    final double d7 = d5 * zMatrix.getEntry(i, 0) + d6 * zMatrix.getEntry(i, j);
+                    zMatrix.setEntry(i, j, d5 * zMatrix.getEntry(i, j) - d6 * zMatrix.getEntry(i, 0));
+                    zMatrix.setEntry(i, 0, d7);
                 }
             }
-            zmat.setEntry(knew, j, ZERO);
+            zMatrix.setEntry(knew, j, ZERO);
         }
 
         // Put the first NPT components of the KNEW-th column of HLAG into W,
         // and calculate the parameters of the updating formula.
 
         for (int i = 0; i < npt; i++) {
-            work.setEntry(i, zmat.getEntry(knew, 0) * zmat.getEntry(i, 0));
+            work.setEntry(i, zMatrix.getEntry(knew, 0) * zMatrix.getEntry(i, 0));
         }
         final double alpha = work.getEntry(knew);
-        final double tau = vlag.getEntry(knew);
-        vlag.setEntry(knew, vlag.getEntry(knew) - ONE);
+        final double tau = lagrangeValuesAtNewPoint.getEntry(knew);
+        lagrangeValuesAtNewPoint.setEntry(knew, lagrangeValuesAtNewPoint.getEntry(knew) - ONE);
 
         // Complete the updating of ZMAT.
 
         final double sqrtDenom = Math.sqrt(denom);
         final double d1 = tau / sqrtDenom;
-        final double d2 = zmat.getEntry(knew, 0) / sqrtDenom;
+        final double d2 = zMatrix.getEntry(knew, 0) / sqrtDenom;
         for (int i = 0; i < npt; i++) {
-            zmat.setEntry(i, 0,
-                          d1 * zmat.getEntry(i, 0) - d2 * vlag.getEntry(i));
+            zMatrix.setEntry(i, 0,
+                          d1 * zMatrix.getEntry(i, 0) - d2 * lagrangeValuesAtNewPoint.getEntry(i));
         }
 
         // Finally, update the matrix BMAT.
 
         for (int j = 0; j < n; j++) {
             final int jp = npt + j;
-            work.setEntry(jp, bmat.getEntry(knew, j));
-            final double d3 = (alpha * vlag.getEntry(jp) - tau * work.getEntry(jp)) / denom;
-            final double d4 = (-beta * work.getEntry(jp) - tau * vlag.getEntry(jp)) / denom;
+            work.setEntry(jp, bMatrix.getEntry(knew, j));
+            final double d3 = (alpha * lagrangeValuesAtNewPoint.getEntry(jp) - tau * work.getEntry(jp)) / denom;
+            final double d4 = (-beta * work.getEntry(jp) - tau * lagrangeValuesAtNewPoint.getEntry(jp)) / denom;
             for (int i = 0; i <= jp; i++) {
-                bmat.setEntry(i, j,
-                              bmat.getEntry(i, j) + d3 * vlag.getEntry(i) + d4 * work.getEntry(i));
+                bMatrix.setEntry(i, j,
+                              bMatrix.getEntry(i, j) + d3 * lagrangeValuesAtNewPoint.getEntry(i) + d4 * work.getEntry(i));
                 if (i >= npt) {
-                    bmat.setEntry(jp, (i - npt), bmat.getEntry(i, j));
+                    bMatrix.setEntry(jp, (i - npt), bMatrix.getEntry(i, j));
                 }
             }
         }
@@ -2421,7 +2359,7 @@ public class BOBYQAOptimizer
      * {@link #upperBound} array if no constraints were provided.
      */
     private void setup() {
-        // System.out.println("setup"); // XXX
+        printMethod(); // XXX
 
         double[] init = getStartPoint();
         final int dimension = init.length;
@@ -2474,6 +2412,26 @@ public class BOBYQAOptimizer
         if (minDiff < requiredMinDiff) {
             initialTrustRegionRadius = minDiff / 3.0;
         }
+
+        // Initialize the data structures used by the "bobyqa" method.
+        bMatrix = new Array2DRowRealMatrix(dimension + numberOfInterpolationPoints,
+                                           dimension);

[... 43 lines stripped ...]


Mime
View raw message