An implementation of Ruhe & Wedin's Algorithm II.

A non-linear generalization of Ruhe & Wedin's algorithm
for separable non-linear least squares problem. It is implemented
as coordinate descent on an independent subset of the parameter
blocks at the end of every successful Newton step. The resulting
algorithm has much improved convergence at the cost of some
execution time.

Change-Id: I8fdc5edbd0ba1e702c9658b98041b2c2ae705402
diff --git a/internal/ceres/trust_region_minimizer.cc b/internal/ceres/trust_region_minimizer.cc
index dd49f9e..b838b19 100644
--- a/internal/ceres/trust_region_minimizer.cc
+++ b/internal/ceres/trust_region_minimizer.cc
@@ -45,6 +45,7 @@
 #include "ceres/internal/scoped_ptr.h"
 #include "ceres/linear_least_squares_problems.h"
 #include "ceres/sparse_matrix.h"
+#include "ceres/stringprintf.h"
 #include "ceres/trust_region_strategy.h"
 #include "ceres/types.h"
 #include "ceres/wall_time.h"
@@ -102,8 +103,6 @@
   //
   // Both of these indicate that this is the wrong place for this
   // code, and going forward this should needs fixing/refactoring.
-  LOG(WARNING) << "Dumping linear least squares problem to disk is "
-               << "currently broken.";
   return true;
 }
 
@@ -142,7 +141,6 @@
   iteration_summary.iteration = 0;
   iteration_summary.step_is_valid = false;
   iteration_summary.step_is_successful = false;
-  iteration_summary.cost = summary->initial_cost;
   iteration_summary.cost_change = 0.0;
   iteration_summary.gradient_max_norm = 0.0;
   iteration_summary.step_norm = 0.0;
@@ -162,6 +160,8 @@
     return;
   }
 
+  iteration_summary.cost = cost + summary->fixed_cost;
+
   int num_consecutive_nonmonotonic_steps = 0;
   double minimum_cost = cost;
   double reference_cost = cost;
@@ -294,10 +294,12 @@
       if (++num_consecutive_invalid_steps >=
           options_.max_num_consecutive_invalid_steps) {
         summary->termination_type = NUMERICAL_FAILURE;
-        LOG(WARNING) << "Terminating. Number of successive invalid steps more "
-                     << "than "
-                     << "Solver::Options::max_num_consecutive_invalid_steps: "
-                     << options_.max_num_consecutive_invalid_steps;
+        summary->error = StringPrintf(
+            "Terminating. Number of successive invalid steps more "
+            "than Solver::Options::max_num_consecutive_invalid_steps: %d",
+            options_.max_num_consecutive_invalid_steps);
+
+        LOG(WARNING) << summary->error;
         return;
       }
 
@@ -306,7 +308,7 @@
       // as an unsuccessful iteration. Since the various callbacks are
       // still executed, we are going to fill the iteration summary
       // with data that assumes a step of length zero and no progress.
-      iteration_summary.cost = cost;
+      iteration_summary.cost = cost + summary->fixed_cost;
       iteration_summary.cost_change = 0.0;
       iteration_summary.gradient_max_norm =
           summary->iterations.back().gradient_max_norm;
@@ -319,7 +321,33 @@
 
       // Undo the Jacobian column scaling.
       delta = (trust_region_step.array() * scale.array()).matrix();
-      iteration_summary.step_norm = delta.norm();
+      if (!evaluator->Plus(x.data(), delta.data(), x_plus_delta.data())) {
+        summary->termination_type = NUMERICAL_FAILURE;
+        summary->error =
+            "Terminating. Failed to compute Plus(x, delta, x_plus_delta).";
+
+        LOG(WARNING) << summary->error;
+        return;
+      }
+
+      // Try this step.
+      double new_cost = numeric_limits<double>::max();
+      if (!evaluator->Evaluate(x_plus_delta.data(),
+                               &new_cost,
+                               NULL, NULL, NULL)) {
+        // If the evaluation of the new cost fails, treat it as a step
+        // with high cost.
+        LOG(WARNING) << "Step failed to evaluate. "
+                     << "Treating it as step with infinite cost";
+        new_cost = numeric_limits<double>::max();
+      } else if (new_cost < cost && options.inner_iteration_minimizer != NULL) {
+        Solver::Summary inner_iteration_summary;
+        options.inner_iteration_minimizer->Minimize(options,
+                                                    x_plus_delta.data(),
+                                                    &inner_iteration_summary);
+      }
+
+      iteration_summary.step_norm = (x - x_plus_delta).norm();
 
       // Convergence based on parameter_tolerance.
       const double step_size_tolerance =  options_.parameter_tolerance *
@@ -334,25 +362,6 @@
         return;
       }
 
-      if (!evaluator->Plus(x.data(), delta.data(), x_plus_delta.data())) {
-        summary->termination_type = NUMERICAL_FAILURE;
-        LOG(WARNING) << "Terminating. Failed to compute "
-                     << "Plus(x, delta, x_plus_delta).";
-        return;
-      }
-
-      // Try this step.
-      double new_cost;
-      if (!evaluator->Evaluate(x_plus_delta.data(),
-                               &new_cost,
-                               NULL, NULL, NULL)) {
-        // If the evaluation of the new cost fails, treat it as a step
-        // with high cost.
-        LOG(WARNING) << "Step failed to evaluate. "
-                     << "Treating it as step with infinite cost";
-        new_cost = numeric_limits<double>::max();
-      }
-
       VLOG(2) << "old cost: " << cost << " new cost: " << new_cost;
       iteration_summary.cost_change =  cost - new_cost;
       const double absolute_function_tolerance =
@@ -420,7 +429,8 @@
                                NULL,
                                jacobian)) {
         summary->termination_type = NUMERICAL_FAILURE;
-        LOG(WARNING) << "Terminating: Residual and Jacobian evaluation failed.";
+        summary->error = "Terminating: Residual and Jacobian evaluation failed.";
+        LOG(WARNING) << summary->error;
         return;
       }