Changes to how gradient based convergence is diagnosed.

The original implementation for computing the norm of the gradient was

gradient_norm = norm(gradient)

when the gradient vector lies in the same space as the parameter
vector, this value is meaningful. When there is a local parameterization
involved, interpreting this value and diagnosing convergence using it
is hard.

Further, this expression does not respect the bounds constraints
on the parmeters. Measuring the norm of the gradient only makes
sense when the optimization being performed is unconstrained.

A better solution, used by LANCELOT is the expression
gradient_norm = norm(x - P(x - gradient))

Here, P is the projection operator onto the bounds constraints.
x - gradient is computed by computing Plus(x, -gradient), thus the
actual expression becomes

gradient_norm = norm(x - P(Plus(x, -gradient)));

Which in the case where there are no bounds constraints, and there
are no local parameterizations, reduces to the usual Euclidean
expression from above, since

Plus(x, -gradient) = x - gradient

and P(x - gradient) = x - gradient.

This change implements this change. Further, the convergence
test using the gradient tolerance now uses an absolute measure
rather than a relative measure. This is a forward looking change
as we start implementing the Augmented Lagrangian solver.

Last but not the least, various "Terminating: Foo" messages
have been changed so that "Terminating: " is logged but is
not part of the Solver::Summary::message string as it is
pointless.

Change-Id: I943146f71a1da47c8c7592986039b4112781b99b
diff --git a/internal/ceres/line_search_minimizer.cc b/internal/ceres/line_search_minimizer.cc
index 339d275..977fd6c 100644
--- a/internal/ceres/line_search_minimizer.cc
+++ b/internal/ceres/line_search_minimizer.cc
@@ -74,18 +74,29 @@
 // this more carefully.
 bool Evaluate(Evaluator* evaluator,
               const Vector& x,
-              LineSearchMinimizer::State* state) {
-  const bool status = evaluator->Evaluate(x.data(),
-                                          &(state->cost),
-                                          NULL,
-                                          state->gradient.data(),
-                                          NULL);
-  if (status) {
-    state->gradient_squared_norm = state->gradient.squaredNorm();
-    state->gradient_max_norm = state->gradient.lpNorm<Eigen::Infinity>();
+              LineSearchMinimizer::State* state,
+              string* message) {
+  if (!evaluator->Evaluate(x.data(),
+                           &(state->cost),
+                           NULL,
+                           state->gradient.data(),
+                           NULL)) {
+    *message = "Gradient evaluation failed.";
+    return false;
   }
 
-  return status;
+  Vector negative_gradient = -state->gradient;
+  Vector projected_gradient_step(x.size());
+  if (!evaluator->Plus(x.data(),
+                       negative_gradient.data(),
+                       projected_gradient_step.data())) {
+    *message = "projected_gradient_step = Plus(x, -gradient) failed.";
+    return false;
+  }
+
+  state->gradient_squared_norm = (x - projected_gradient_step).squaredNorm();
+  state->gradient_max_norm = (x - projected_gradient_step).lpNorm<Eigen::Infinity>();
+  return true;
 }
 
 }  // namespace
@@ -125,10 +136,9 @@
   iteration_summary.step_solver_time_in_seconds = 0;
 
   // Do initial cost and Jacobian evaluation.
-  if (!Evaluate(evaluator, x, &current_state)) {
-    summary->message = "Terminating: Cost and gradient evaluation failed.";
+  if (!Evaluate(evaluator, x, &current_state, &summary->message)) {
     summary->termination_type = FAILURE;
-    LOG_IF(WARNING, is_not_silent) << summary->message;
+    LOG_IF(WARNING, is_not_silent) << "Terminating: " << summary->message;
     return;
   }
 
@@ -138,22 +148,13 @@
   iteration_summary.gradient_max_norm = current_state.gradient_max_norm;
   iteration_summary.gradient_norm = sqrt(current_state.gradient_squared_norm);
 
-  // The initial gradient max_norm is bounded from below so that we do
-  // not divide by zero.
-  const double initial_gradient_max_norm =
-      max(iteration_summary.gradient_max_norm, kEpsilon);
-  const double absolute_gradient_tolerance =
-      options.gradient_tolerance * initial_gradient_max_norm;
-
-  if (iteration_summary.gradient_max_norm <= absolute_gradient_tolerance) {
-    summary->message =
-        StringPrintf("Terminating: Gradient tolerance reached. "
-                     "Relative gradient max norm: %e <= %e",
-                     iteration_summary.gradient_max_norm /
-                     initial_gradient_max_norm,
-                     options.gradient_tolerance);
+  if (iteration_summary.gradient_max_norm <= options.gradient_tolerance) {
+    summary->message = StringPrintf("Gradient tolerance reached. "
+                                    "Gradient max norm: %e <= %e",
+                                    iteration_summary.gradient_max_norm,
+                                    options.gradient_tolerance);
     summary->termination_type = CONVERGENCE;
-    VLOG_IF(1, is_not_silent) << summary->message;
+    VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
     return;
   }
 
@@ -210,23 +211,23 @@
 
   while (true) {
     if (!RunCallbacks(options.callbacks, iteration_summary, summary)) {
-      return;
+      break;
     }
 
     iteration_start_time = WallTimeInSeconds();
     if (iteration_summary.iteration >= options.max_num_iterations) {
-      summary->message = "Terminating: Maximum number of iterations reached.";
+      summary->message = "Maximum number of iterations reached.";
       summary->termination_type = NO_CONVERGENCE;
-      VLOG_IF(1, is_not_silent) << summary->message;
+      VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
       break;
     }
 
     const double total_solver_time = iteration_start_time - start_time +
         summary->preprocessor_time_in_seconds;
     if (total_solver_time >= options.max_solver_time_in_seconds) {
-      summary->message = "Terminating: Maximum solver time reached.";
+      summary->message = "Maximum solver time reached.";
       summary->termination_type = NO_CONVERGENCE;
-      VLOG_IF(1, is_not_silent) << summary->message;
+      VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
       break;
     }
 
@@ -252,11 +253,11 @@
       // have already reached our specified maximum number of restarts,
       // terminate optimization.
       summary->message =
-          StringPrintf("Terminating: Line search direction failure: specified "
+          StringPrintf("Line search direction failure: specified "
                        "max_num_line_search_direction_restarts: %d reached.",
                        options.max_num_line_search_direction_restarts);
       summary->termination_type = FAILURE;
-      LOG_IF(WARNING, is_not_silent) << summary->message;
+      LOG_IF(WARNING, is_not_silent) << "Terminating: " << summary->message;
       break;
     } else if (!line_search_status) {
       // Restart line search direction with gradient descent on first iteration
@@ -345,8 +346,12 @@
     if (!evaluator->Plus(x.data(), delta.data(), x_plus_delta.data())) {
       LOG_IF(WARNING, is_not_silent)
           << "x_plus_delta = Plus(x, delta) failed. ";
-    } else if (!Evaluate(evaluator, x_plus_delta, &current_state)) {
-      LOG_IF(WARNING, is_not_silent) << "Step failed to evaluate. ";
+    } else if (!Evaluate(evaluator,
+                         x_plus_delta,
+                         &current_state,
+                         &summary->message)) {
+      LOG_IF(WARNING, is_not_silent) << "Step failed to evaluate. "
+                                     << summary->message;
     } else {
       x = x_plus_delta;
     }
@@ -354,15 +359,13 @@
     iteration_summary.gradient_max_norm = current_state.gradient_max_norm;
     iteration_summary.gradient_norm = sqrt(current_state.gradient_squared_norm);
 
-    if (iteration_summary.gradient_max_norm <= absolute_gradient_tolerance) {
-      summary->message =
-          StringPrintf("Terminating: Gradient tolerance reached. "
-                       "Relative gradient max norm: %e <= %e. ",
-                       (iteration_summary.gradient_max_norm /
-                        initial_gradient_max_norm),
-                       options.gradient_tolerance);
+    if (iteration_summary.gradient_max_norm <= options.gradient_tolerance) {
+      summary->message = StringPrintf("Gradient tolerance reached. "
+                                      "Gradient max norm: %e <= %e",
+                                      iteration_summary.gradient_max_norm,
+                                      options.gradient_tolerance);
       summary->termination_type = CONVERGENCE;
-      VLOG_IF(1, is_not_silent) << summary->message;
+      VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
       break;
     }
 
@@ -371,14 +374,14 @@
         options.function_tolerance * previous_state.cost;
     if (fabs(iteration_summary.cost_change) < absolute_function_tolerance) {
       summary->message =
-          StringPrintf("Terminating. Function tolerance reached. "
+          StringPrintf("Function tolerance reached. "
                        "|cost_change|/cost: %e <= %e",
                        fabs(iteration_summary.cost_change) /
                        previous_state.cost,
                        options.function_tolerance);
       summary->termination_type = CONVERGENCE;
-      VLOG_IF(1, is_not_silent) << summary->message;
-      return;
+      VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
+      break;
     }
 
     iteration_summary.cost = current_state.cost + summary->fixed_cost;
diff --git a/internal/ceres/solver.cc b/internal/ceres/solver.cc
index da87de1..36c58b3 100644
--- a/internal/ceres/solver.cc
+++ b/internal/ceres/solver.cc
@@ -344,8 +344,8 @@
   StringAppendF(&report, "Total               %25.3f\n\n",
                 total_time_in_seconds);
 
-  StringAppendF(&report, "Termination:        %25s\n",
-                TerminationTypeToString(termination_type));
+  StringAppendF(&report, "Termination:        %25s (%s)\n",
+                TerminationTypeToString(termination_type), message.c_str());
   return report;
 };
 
diff --git a/internal/ceres/trust_region_minimizer.cc b/internal/ceres/trust_region_minimizer.cc
index f8ecac2..51eb75f 100644
--- a/internal/ceres/trust_region_minimizer.cc
+++ b/internal/ceres/trust_region_minimizer.cc
@@ -152,6 +152,8 @@
   Vector gradient(num_effective_parameters);
   Vector model_residuals(num_residuals);
   Vector scale(num_effective_parameters);
+  Vector negative_gradient(num_effective_parameters);
+  Vector projected_gradient_step(num_parameters);
 
   IterationSummary iteration_summary;
   iteration_summary.iteration = 0;
@@ -197,36 +199,39 @@
     return;
   }
 
-  int num_consecutive_nonmonotonic_steps = 0;
-  double minimum_cost = cost;
-  double reference_cost = cost;
-  double accumulated_reference_model_cost_change = 0.0;
-  double candidate_cost = cost;
-  double accumulated_candidate_model_cost_change = 0.0;
+  negative_gradient = -gradient;
+  if (!evaluator->Plus(x.data(),
+                       negative_gradient.data(),
+                       projected_gradient_step.data())) {
+    summary->message = "Unable to compute gradient step.";
+    summary->termination_type = FAILURE;
+    LOG(ERROR) << "Terminating: " << summary->message;
+    return;
+  }
 
   summary->initial_cost = cost + summary->fixed_cost;
   iteration_summary.cost = cost + summary->fixed_cost;
-  iteration_summary.gradient_max_norm = gradient.lpNorm<Eigen::Infinity>();
-  iteration_summary.gradient_norm = gradient.norm();
+  iteration_summary.gradient_max_norm =
+    (x - projected_gradient_step).lpNorm<Eigen::Infinity>();
+  iteration_summary.gradient_norm = (x - projected_gradient_step).norm();
 
-  // The initial gradient max_norm is bounded from below so that we do
-  // not divide by zero.
-  const double initial_gradient_max_norm =
-      max(iteration_summary.gradient_max_norm, kEpsilon);
-  const double absolute_gradient_tolerance =
-      options_.gradient_tolerance * initial_gradient_max_norm;
-
-  if (iteration_summary.gradient_max_norm <= absolute_gradient_tolerance) {
-    summary->message = StringPrintf("Terminating: Gradient tolerance reached. "
-                                  "Relative gradient max norm: %e <= %e",
-                                  (iteration_summary.gradient_max_norm /
-                                   initial_gradient_max_norm),
-                                  options_.gradient_tolerance);
+  if (iteration_summary.gradient_max_norm <= options.gradient_tolerance) {
+    summary->message = StringPrintf("Gradient tolerance reached. "
+                                    "Gradient max norm: %e <= %e",
+                                    iteration_summary.gradient_max_norm,
+                                    options_.gradient_tolerance);
     summary->termination_type = CONVERGENCE;
-    VLOG_IF(1, is_not_silent) << summary->message;
+    VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
     return;
   }
 
+  if (options_.jacobi_scaling) {
+    EstimateScale(*jacobian, scale.data());
+    jacobian->ScaleColumns(scale.data());
+  } else {
+    scale.setOnes();
+  }
+
   iteration_summary.iteration_time_in_seconds =
       WallTimeInSeconds() - iteration_start_time;
   iteration_summary.cumulative_time_in_seconds =
@@ -234,13 +239,12 @@
       + summary->preprocessor_time_in_seconds;
   summary->iterations.push_back(iteration_summary);
 
-  if (options_.jacobi_scaling) {
-    EstimateScale(*jacobian, scale.data());
-    jacobian->ScaleColumns(scale.data());
-  } else {
-    scale.setOnes();
-  }
-
+  int num_consecutive_nonmonotonic_steps = 0;
+  double minimum_cost = cost;
+  double reference_cost = cost;
+  double accumulated_reference_model_cost_change = 0.0;
+  double candidate_cost = cost;
+  double accumulated_candidate_model_cost_change = 0.0;
   int num_consecutive_invalid_steps = 0;
   bool inner_iterations_are_enabled = options.inner_iteration_minimizer != NULL;
   while (true) {
@@ -251,18 +255,18 @@
 
     iteration_start_time = WallTimeInSeconds();
     if (iteration_summary.iteration >= options_.max_num_iterations) {
-      summary->message = "Terminating: Maximum number of iterations reached.";
+      summary->message = "Maximum number of iterations reached.";
       summary->termination_type = NO_CONVERGENCE;
-      VLOG_IF(1, is_not_silent) << summary->message;
+      VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
       return;
     }
 
     const double total_solver_time = iteration_start_time - start_time +
         summary->preprocessor_time_in_seconds;
     if (total_solver_time >= options_.max_solver_time_in_seconds) {
-      summary->message = "Terminating: Maximum solver time reached.";
+      summary->message = "Maximum solver time reached.";
       summary->termination_type = NO_CONVERGENCE;
-      VLOG_IF(1, is_not_silent) << summary->message;
+      VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
       return;
     }
 
@@ -292,10 +296,10 @@
 
     if (strategy_summary.termination_type == LINEAR_SOLVER_FATAL_ERROR) {
       summary->message =
-          "Terminating. Linear solver failed due to unrecoverable "
+          "Linear solver failed due to unrecoverable "
           "non-numeric causes. Please see the error log for clues. ";
       summary->termination_type = FAILURE;
-      LOG_IF(WARNING, is_not_silent) << summary->message;
+      LOG_IF(WARNING, is_not_silent) << "Terminating: " << summary->message;
       return;
     }
 
@@ -339,11 +343,11 @@
       if (++num_consecutive_invalid_steps >=
           options_.max_num_consecutive_invalid_steps) {
         summary->message = StringPrintf(
-            "Terminating. Number of successive invalid steps more "
+            "Number of successive invalid steps more "
             "than Solver::Options::max_num_consecutive_invalid_steps: %d",
             options_.max_num_consecutive_invalid_steps);
         summary->termination_type = FAILURE;
-        LOG_IF(WARNING, is_not_silent) << summary->message;
+        LOG_IF(WARNING, is_not_silent) << "Terminating: " << summary->message;
         return;
       }
 
@@ -448,13 +452,13 @@
           (x_norm + options_.parameter_tolerance);
       if (iteration_summary.step_norm <= step_size_tolerance) {
         summary->message =
-            StringPrintf("Terminating. Parameter tolerance reached. "
-                         "relative step_norm: %e <= %e.",
+            StringPrintf("Parameter tolerance reached. "
+                         "Relative step_norm: %e <= %e.",
                          (iteration_summary.step_norm /
                           (x_norm + options_.parameter_tolerance)),
                          options_.parameter_tolerance);
         summary->termination_type = CONVERGENCE;
-        VLOG_IF(1, is_not_silent) << summary->message;
+        VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
         return;
       }
 
@@ -463,12 +467,12 @@
           options_.function_tolerance * cost;
       if (fabs(iteration_summary.cost_change) < absolute_function_tolerance) {
         summary->message =
-            StringPrintf("Terminating. Function tolerance reached. "
+            StringPrintf("Function tolerance reached. "
                          "|cost_change|/cost: %e <= %e",
                          fabs(iteration_summary.cost_change) / cost,
                          options_.function_tolerance);
         summary->termination_type = CONVERGENCE;
-        VLOG_IF(1, is_not_silent) << summary->message;
+        VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
         return;
       }
 
@@ -565,25 +569,34 @@
                                residuals.data(),
                                gradient.data(),
                                jacobian)) {
-        summary->message =
-            "Terminating: Residual and Jacobian evaluation failed.";
+        summary->message = "Residual and Jacobian evaluation failed.";
         summary->termination_type = FAILURE;
-        LOG_IF(WARNING, is_not_silent) << summary->message;
+        LOG_IF(WARNING, is_not_silent) << "Terminating: " << summary->message;
         return;
       }
 
-      iteration_summary.gradient_max_norm = gradient.lpNorm<Eigen::Infinity>();
-      iteration_summary.gradient_norm = gradient.norm();
-
-      if (iteration_summary.gradient_max_norm <= absolute_gradient_tolerance) {
+      negative_gradient = -gradient;
+      if (!evaluator->Plus(x.data(),
+                           negative_gradient.data(),
+                           projected_gradient_step.data())) {
         summary->message =
-            StringPrintf("Terminating: Gradient tolerance reached. "
-                         "Relative gradient max norm: %e <= %e",
-                         (iteration_summary.gradient_max_norm /
-                          initial_gradient_max_norm),
-                         options_.gradient_tolerance);
+            "projected_gradient_step = Plus(x, -gradient) failed.";
+        summary->termination_type = FAILURE;
+        LOG(ERROR) << "Terminating: " << summary->message;
+        return;
+      }
+
+      iteration_summary.gradient_max_norm =
+        (x - projected_gradient_step).lpNorm<Eigen::Infinity>();
+      iteration_summary.gradient_norm = (x - projected_gradient_step).norm();
+
+      if (iteration_summary.gradient_max_norm <= options.gradient_tolerance) {
+        summary->message = StringPrintf("Gradient tolerance reached. "
+                                        "Gradient max norm: %e <= %e",
+                                        iteration_summary.gradient_max_norm,
+                                        options_.gradient_tolerance);
         summary->termination_type = CONVERGENCE;
-        VLOG_IF(1, is_not_silent) << summary->message;
+        VLOG_IF(1, is_not_silent) << "Terminating: " << summary->message;
         return;
       }