internal/ceres/numeric_diff_cost_function_test.cc

   1 // Ceres Solver - A fast non-linear least squares minimizer
   2 // Copyright 2015 Google Inc. All rights reserved.
   3 // http://ceres-solver.org/
   4 //
   5 // Redistribution and use in source and binary forms, with or without
   6 // modification, are permitted provided that the following conditions are met:
   7 //
   8 // * Redistributions of source code must retain the above copyright notice,
   9 //   this list of conditions and the following disclaimer.
  10 // * Redistributions in binary form must reproduce the above copyright notice,
  11 //   this list of conditions and the following disclaimer in the documentation
  12 //   and/or other materials provided with the distribution.
  13 // * Neither the name of Google Inc. nor the names of its contributors may be
  14 //   used to endorse or promote products derived from this software without
  15 //   specific prior written permission.
  16 //
  17 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
  18 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  19 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  20 // ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
  21 // LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  22 // CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  23 // SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  24 // INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  25 // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  26 // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  27 // POSSIBILITY OF SUCH DAMAGE.
  28 //
  29 // Author: keir@google.com (Keir Mierle)
  30 //         tbennun@gmail.com (Tal Ben-Nun)
  31
  32 #include "ceres/numeric_diff_cost_function.h"
  33
  34 #include <algorithm>
  35 #include <cmath>
  36 #include <string>
  37 #include <vector>
  38 #include "ceres/internal/macros.h"
  39 #include "ceres/internal/scoped_ptr.h"
  40 #include "ceres/array_utils.h"
  41 #include "ceres/numeric_diff_test_utils.h"
  42 #include "ceres/test_util.h"
  43 #include "ceres/types.h"
  44 #include "glog/logging.h"
  45 #include "gtest/gtest.h"
  46
  47 namespace ceres {
  48 namespace internal {
  49
  50 TEST(NumericDiffCostFunction, EasyCaseFunctorCentralDifferences) {
  51   internal::scoped_ptr<CostFunction> cost_function;
  52   cost_function.reset(
  53       new NumericDiffCostFunction<EasyFunctor,
  54                                   CENTRAL,
  55                                   3,  /* number of residuals */
  56                                   5,  /* size of x1 */
  57                                   5   /* size of x2 */>(
  58           new EasyFunctor));
  59   EasyFunctor functor;
  60   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, CENTRAL);
  61 }
  62
  63 TEST(NumericDiffCostFunction, EasyCaseFunctorForwardDifferences) {
  64   internal::scoped_ptr<CostFunction> cost_function;
  65   cost_function.reset(
  66       new NumericDiffCostFunction<EasyFunctor,
  67                                   FORWARD,
  68                                   3,  /* number of residuals */
  69                                   5,  /* size of x1 */
  70                                   5   /* size of x2 */>(
  71           new EasyFunctor));
  72   EasyFunctor functor;
  73   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, FORWARD);
  74 }
  75
  76 TEST(NumericDiffCostFunction, EasyCaseFunctorRidders) {
  77   internal::scoped_ptr<CostFunction> cost_function;
  78   cost_function.reset(
  79       new NumericDiffCostFunction<EasyFunctor,
  80                                   RIDDERS,
  81                                   3,  /* number of residuals */
  82                                   5,  /* size of x1 */
  83                                   5   /* size of x2 */>(
  84           new EasyFunctor));
  85   EasyFunctor functor;
  86   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, RIDDERS);
  87 }
  88
  89 TEST(NumericDiffCostFunction, EasyCaseCostFunctionCentralDifferences) {
  90   internal::scoped_ptr<CostFunction> cost_function;
  91   cost_function.reset(
  92       new NumericDiffCostFunction<EasyCostFunction,
  93                                   CENTRAL,
  94                                   3,  /* number of residuals */
  95                                   5,  /* size of x1 */
  96                                   5   /* size of x2 */>(
  97           new EasyCostFunction, TAKE_OWNERSHIP));
  98   EasyFunctor functor;
  99   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, CENTRAL);
 100 }
 101
 102 TEST(NumericDiffCostFunction, EasyCaseCostFunctionForwardDifferences) {
 103   internal::scoped_ptr<CostFunction> cost_function;
 104   cost_function.reset(
 105       new NumericDiffCostFunction<EasyCostFunction,
 106                                   FORWARD,
 107                                   3,  /* number of residuals */
 108                                   5,  /* size of x1 */
 109                                   5   /* size of x2 */>(
 110           new EasyCostFunction, TAKE_OWNERSHIP));
 111   EasyFunctor functor;
 112   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, FORWARD);
 113 }
 114
 115 TEST(NumericDiffCostFunction, EasyCaseCostFunctionRidders) {
 116   internal::scoped_ptr<CostFunction> cost_function;
 117   cost_function.reset(
 118       new NumericDiffCostFunction<EasyCostFunction,
 119                                   RIDDERS,
 120                                   3,  /* number of residuals */
 121                                   5,  /* size of x1 */
 122                                   5   /* size of x2 */>(
 123           new EasyCostFunction, TAKE_OWNERSHIP));
 124   EasyFunctor functor;
 125   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, RIDDERS);
 126 }
 127
 128 TEST(NumericDiffCostFunction,
 129      TranscendentalCaseFunctorCentralDifferences) {
 130   internal::scoped_ptr<CostFunction> cost_function;
 131   cost_function.reset(
 132       new NumericDiffCostFunction<TranscendentalFunctor,
 133                                   CENTRAL,
 134                                   2,  /* number of residuals */
 135                                   5,  /* size of x1 */
 136                                   5   /* size of x2 */>(
 137           new TranscendentalFunctor));
 138   TranscendentalFunctor functor;
 139   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, CENTRAL);
 140 }
 141
 142 TEST(NumericDiffCostFunction,
 143      TranscendentalCaseFunctorForwardDifferences) {
 144   internal::scoped_ptr<CostFunction> cost_function;
 145   cost_function.reset(
 146       new NumericDiffCostFunction<TranscendentalFunctor,
 147                                   FORWARD,
 148                                   2,  /* number of residuals */
 149                                   5,  /* size of x1 */
 150                                   5   /* size of x2 */>(
 151           new TranscendentalFunctor));
 152   TranscendentalFunctor functor;
 153   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, FORWARD);
 154 }
 155
 156 TEST(NumericDiffCostFunction,
 157      TranscendentalCaseFunctorRidders) {
 158   NumericDiffOptions options;
 159
 160   // Using a smaller initial step size to overcome oscillatory function
 161   // behavior.
 162   options.ridders_relative_initial_step_size = 1e-3;
 163
 164   internal::scoped_ptr<CostFunction> cost_function;
 165   cost_function.reset(
 166       new NumericDiffCostFunction<TranscendentalFunctor,
 167                                   RIDDERS,
 168                                   2,  /* number of residuals */
 169                                   5,  /* size of x1 */
 170                                   5   /* size of x2 */>(
 171           new TranscendentalFunctor, TAKE_OWNERSHIP, 2, options));
 172   TranscendentalFunctor functor;
 173   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, RIDDERS);
 174 }
 175
 176 TEST(NumericDiffCostFunction,
 177      TranscendentalCaseCostFunctionCentralDifferences) {
 178   internal::scoped_ptr<CostFunction> cost_function;
 179   cost_function.reset(
 180       new NumericDiffCostFunction<TranscendentalCostFunction,
 181                                   CENTRAL,
 182                                   2,  /* number of residuals */
 183                                   5,  /* size of x1 */
 184                                   5   /* size of x2 */>(
 185           new TranscendentalCostFunction, TAKE_OWNERSHIP));
 186   TranscendentalFunctor functor;
 187   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, CENTRAL);
 188 }
 189
 190 TEST(NumericDiffCostFunction,
 191      TranscendentalCaseCostFunctionForwardDifferences) {
 192   internal::scoped_ptr<CostFunction> cost_function;
 193   cost_function.reset(
 194       new NumericDiffCostFunction<TranscendentalCostFunction,
 195                                   FORWARD,
 196                                   2,  /* number of residuals */
 197                                   5,  /* size of x1 */
 198                                   5   /* size of x2 */>(
 199           new TranscendentalCostFunction, TAKE_OWNERSHIP));
 200   TranscendentalFunctor functor;
 201   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, FORWARD);
 202 }
 203
 204 TEST(NumericDiffCostFunction,
 205      TranscendentalCaseCostFunctionRidders) {
 206   NumericDiffOptions options;
 207
 208   // Using a smaller initial step size to overcome oscillatory function
 209   // behavior.
 210   options.ridders_relative_initial_step_size = 1e-3;
 211
 212   internal::scoped_ptr<CostFunction> cost_function;
 213   cost_function.reset(
 214       new NumericDiffCostFunction<TranscendentalCostFunction,
 215                                   RIDDERS,
 216                                   2,  /* number of residuals */
 217                                   5,  /* size of x1 */
 218                                   5   /* size of x2 */>(
 219           new TranscendentalCostFunction, TAKE_OWNERSHIP, 2, options));
 220   TranscendentalFunctor functor;
 221   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, RIDDERS);
 222 }
 223
 224 template<int num_rows, int num_cols>
 225 class SizeTestingCostFunction : public SizedCostFunction<num_rows, num_cols> {
 226  public:
 227   virtual bool Evaluate(double const* const* parameters,
 228                         double* residuals,
 229                         double** jacobians) const {
 230     return true;
 231   }
 232 };
 233
 234 // As described in
 235 // http://forum.kde.org/viewtopic.php?f=74&t=98536#p210774
 236 // Eigen3 has restrictions on the Row/Column major storage of vectors,
 237 // depending on their dimensions. This test ensures that the correct
 238 // templates are instantiated for various shapes of the Jacobian
 239 // matrix.
 240 TEST(NumericDiffCostFunction, EigenRowMajorColMajorTest) {
 241   scoped_ptr<CostFunction> cost_function;
 242   cost_function.reset(
 243       new NumericDiffCostFunction<SizeTestingCostFunction<1,1>,  CENTRAL, 1, 1>(
 244           new SizeTestingCostFunction<1,1>, ceres::TAKE_OWNERSHIP));
 245
 246   cost_function.reset(
 247       new NumericDiffCostFunction<SizeTestingCostFunction<2,1>,  CENTRAL, 2, 1>(
 248           new SizeTestingCostFunction<2,1>, ceres::TAKE_OWNERSHIP));
 249
 250   cost_function.reset(
 251       new NumericDiffCostFunction<SizeTestingCostFunction<1,2>,  CENTRAL, 1, 2>(
 252           new SizeTestingCostFunction<1,2>, ceres::TAKE_OWNERSHIP));
 253
 254   cost_function.reset(
 255       new NumericDiffCostFunction<SizeTestingCostFunction<2,2>,  CENTRAL, 2, 2>(
 256           new SizeTestingCostFunction<2,2>, ceres::TAKE_OWNERSHIP));
 257
 258   cost_function.reset(
 259       new NumericDiffCostFunction<EasyFunctor, CENTRAL, ceres::DYNAMIC, 1, 1>(
 260           new EasyFunctor, TAKE_OWNERSHIP, 1));
 261
 262   cost_function.reset(
 263       new NumericDiffCostFunction<EasyFunctor, CENTRAL, ceres::DYNAMIC, 1, 1>(
 264           new EasyFunctor, TAKE_OWNERSHIP, 2));
 265
 266   cost_function.reset(
 267       new NumericDiffCostFunction<EasyFunctor, CENTRAL, ceres::DYNAMIC, 1, 2>(
 268           new EasyFunctor, TAKE_OWNERSHIP, 1));
 269
 270   cost_function.reset(
 271       new NumericDiffCostFunction<EasyFunctor, CENTRAL, ceres::DYNAMIC, 1, 2>(
 272           new EasyFunctor, TAKE_OWNERSHIP, 2));
 273
 274   cost_function.reset(
 275       new NumericDiffCostFunction<EasyFunctor, CENTRAL, ceres::DYNAMIC, 2, 1>(
 276           new EasyFunctor, TAKE_OWNERSHIP, 1));
 277
 278   cost_function.reset(
 279       new NumericDiffCostFunction<EasyFunctor, CENTRAL, ceres::DYNAMIC, 2, 1>(
 280           new EasyFunctor, TAKE_OWNERSHIP, 2));
 281 }
 282
 283 TEST(NumericDiffCostFunction,
 284      EasyCaseFunctorCentralDifferencesAndDynamicNumResiduals) {
 285   internal::scoped_ptr<CostFunction> cost_function;
 286   cost_function.reset(
 287       new NumericDiffCostFunction<EasyFunctor,
 288                                   CENTRAL,
 289                                   ceres::DYNAMIC,
 290                                   5,  /* size of x1 */
 291                                   5   /* size of x2 */>(
 292                                       new EasyFunctor, TAKE_OWNERSHIP, 3));
 293   EasyFunctor functor;
 294   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function, CENTRAL);
 295 }
 296
 297 TEST(NumericDiffCostFunction, ExponentialFunctorRidders) {
 298   internal::scoped_ptr<CostFunction> cost_function;
 299   cost_function.reset(
 300       new NumericDiffCostFunction<ExponentialFunctor,
 301                                   RIDDERS,
 302                                   1,  /* number of residuals */
 303                                   1   /* size of x1 */>(
 304              new ExponentialFunctor));
 305   ExponentialFunctor functor;
 306   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function);
 307 }
 308
 309 TEST(NumericDiffCostFunction, ExponentialCostFunctionRidders) {
 310   internal::scoped_ptr<CostFunction> cost_function;
 311   cost_function.reset(
 312       new NumericDiffCostFunction<ExponentialCostFunction,
 313                                   RIDDERS,
 314                                   1,  /* number of residuals */
 315                                   1   /* size of x1 */>(
 316              new ExponentialCostFunction));
 317   ExponentialFunctor functor;
 318   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function);
 319 }
 320
 321 TEST(NumericDiffCostFunction, RandomizedFunctorRidders) {
 322   internal::scoped_ptr<CostFunction> cost_function;
 323   NumericDiffOptions options;
 324   // Larger initial step size is chosen to produce robust results in the
 325   // presence of random noise.
 326   options.ridders_relative_initial_step_size = 10.0;
 327
 328   cost_function.reset(
 329       new NumericDiffCostFunction<RandomizedFunctor,
 330                                   RIDDERS,
 331                                   1,  /* number of residuals */
 332                                   1   /* size of x1 */>(
 333              new RandomizedFunctor(kNoiseFactor, kRandomSeed), TAKE_OWNERSHIP,
 334              1, options));
 335   RandomizedFunctor functor (kNoiseFactor, kRandomSeed);
 336   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function);
 337 }
 338
 339 TEST(NumericDiffCostFunction, RandomizedCostFunctionRidders) {
 340   internal::scoped_ptr<CostFunction> cost_function;
 341   NumericDiffOptions options;
 342   // Larger initial step size is chosen to produce robust results in the
 343   // presence of random noise.
 344   options.ridders_relative_initial_step_size = 10.0;
 345
 346   cost_function.reset(
 347       new NumericDiffCostFunction<RandomizedCostFunction,
 348                                   RIDDERS,
 349                                   1,  /* number of residuals */
 350                                   1   /* size of x1 */>(
 351              new RandomizedCostFunction(kNoiseFactor, kRandomSeed),
 352              TAKE_OWNERSHIP, 1, options));
 353   RandomizedFunctor functor (kNoiseFactor, kRandomSeed);
 354   functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function);
 355 }
 356
 357 struct OnlyFillsOneOutputFunctor {
 358   bool operator()(const double* x, double* output) const {
 359     output[0] = x[0];
 360     return true;
 361   }
 362 };
 363
 364 TEST(NumericDiffCostFunction, PartiallyFilledResidualShouldFailEvaluation) {
 365   double parameter = 1.0;
 366   double jacobian[2];
 367   double residuals[2];
 368   double* parameters[] = {&parameter};
 369   double* jacobians[] = {jacobian};
 370
 371   scoped_ptr<CostFunction> cost_function(
 372       new NumericDiffCostFunction<OnlyFillsOneOutputFunctor, CENTRAL, 2, 1>(
 373           new OnlyFillsOneOutputFunctor));
 374   InvalidateArray(2, jacobian);
 375   InvalidateArray(2, residuals);
 376   EXPECT_TRUE(cost_function->Evaluate(parameters, residuals, jacobians));
 377   EXPECT_FALSE(IsArrayValid(2, residuals));
 378   InvalidateArray(2, residuals);
 379   EXPECT_TRUE(cost_function->Evaluate(parameters, residuals, NULL));
 380   // We are only testing residuals here, because the Jacobians are
 381   // computed using finite differencing from the residuals, so unless
 382   // we introduce a validation step after every evaluation of
 383   // residuals inside NumericDiffCostFunction, there is no way of
 384   // ensuring that the Jacobian array is invalid.
 385   EXPECT_FALSE(IsArrayValid(2, residuals));
 386 }
 387
 388 }  // namespace internal
 389 }  // namespace ceres