← Back to branch summary

~ubuntu-branches/ubuntu/raring/ceres-solver/raring

« back to all changes in this revision

Viewing changes to internal/ceres/numeric_diff_cost_function_test.cc

  • Committer: Package Import Robot
  • Author(s): Koichi Akabe
  • Date: 2012-06-04 07:15:43 UTC
  • Revision ID: package-import@ubuntu.com-20120604071543-zx6uthupvmtqn3k2
Tags: upstream-1.1.1
Import upstream version 1.1.1

Show diffs side-by-side

added added

removed removed

Lines of Context:
internal/ceres/numeric_diff_cost_function_test.cc
 
1
// Ceres Solver - A fast non-linear least squares minimizer
 
2
// Copyright 2010, 2011, 2012 Google Inc. All rights reserved.
 
3
// http://code.google.com/p/ceres-solver/
 
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
 
 
31
#include "ceres/numeric_diff_cost_function.h"
 
32
 
 
33
#include <algorithm>
 
34
#include <cmath>
 
35
#include <string>
 
36
#include <vector>
 
37
#include <glog/logging.h>
 
38
#include "gtest/gtest.h"
 
39
#include "ceres/stringprintf.h"
 
40
#include "ceres/test_util.h"
 
41
#include "ceres/cost_function.h"
 
42
#include "ceres/internal/macros.h"
 
43
#include "ceres/internal/scoped_ptr.h"
 
44
#include "ceres/types.h"
 
45
 
 
46
namespace ceres {
 
47
namespace internal {
 
48
 
 
49
// y1 = x1'x2      -> dy1/dx1 = x2,               dy1/dx2 = x1
 
50
// y2 = (x1'x2)^2  -> dy2/dx1 = 2 * x2 * (x1'x2), dy2/dx2 = 2 * x1 * (x1'x2)
 
51
// y3 = x2'x2      -> dy3/dx1 = 0,                dy3/dx2 = 2 * x2
 
52
class TestCostFunction : public CostFunction {
 
53
 public:
 
54
  TestCostFunction() {
 
55
    set_num_residuals(3);
 
56
    mutable_parameter_block_sizes()->push_back(5);  // x1.
 
57
    mutable_parameter_block_sizes()->push_back(5);  // x2.
 
58
  }
 
59
  virtual bool Evaluate(double const* const* parameters,
 
60
                        double* residuals,
 
61
                        double** jacobians) const {
 
62
    (void) jacobians;  // Ignored.
 
63
 
 
64
    residuals[0] = residuals[1] = residuals[2] = 0;
 
65
    for (int i = 0; i < 5; ++i) {
 
66
      residuals[0] += parameters[0][i] * parameters[1][i];
 
67
      residuals[2] += parameters[1][i] * parameters[1][i];
 
68
    }
 
69
    residuals[1] = residuals[0] * residuals[0];
 
70
    return true;
 
71
  }
 
72
};
 
73
 
 
74
TEST(NumericDiffCostFunction, EasyCase) {
 
75
  // Try both central and forward difference.
 
76
  internal::scoped_ptr<CostFunction> cfs[2];
 
77
  cfs[0].reset(
 
78
      new NumericDiffCostFunction<TestCostFunction,
 
79
                                  CENTRAL,
 
80
                                  3,  /* number of residuals */
 
81
                                  5,  /* size of x1 */
 
82
                                  5   /* size of x2 */>(
 
83
          new TestCostFunction, TAKE_OWNERSHIP));
 
84
 
 
85
  cfs[1].reset(
 
86
      new NumericDiffCostFunction<TestCostFunction,
 
87
                                  FORWARD,
 
88
                                  3,  /* number of residuals */
 
89
                                  5,  /* size of x1 */
 
90
                                  5   /* size of x2 */>(
 
91
          new TestCostFunction, TAKE_OWNERSHIP));
 
92
 
 
93
  for (int c = 0; c < 2; ++c) {
 
94
    CostFunction *cost_function = cfs[c].get();
 
95
 
 
96
    double x1[] = { 1.0, 2.0, 3.0, 4.0, 5.0 };
 
97
    double x2[] = { 9.0, 9.0, 5.0, 5.0, 1.0 };
 
98
    double *parameters[] = { &x1[0], &x2[0] };
 
99
 
 
100
    double dydx1[15];  // 3 x 5, row major.
 
101
    double dydx2[15];  // 3 x 5, row major.
 
102
    double *jacobians[2] = { &dydx1[0], &dydx2[0] };
 
103
 
 
104
    double residuals[3] = {-1e-100, -2e-100, -3e-100 };
 
105
 
 
106
    ASSERT_TRUE(cost_function->Evaluate(&parameters[0],
 
107
                                        &residuals[0],
 
108
                                        &jacobians[0]));
 
109
 
 
110
    EXPECT_EQ(residuals[0], 67);
 
111
    EXPECT_EQ(residuals[1], 4489);
 
112
    EXPECT_EQ(residuals[2], 213);
 
113
 
 
114
    for (int i = 0; i < 5; ++i) {
 
115
      LOG(INFO) << "c = " << c << " i = " << i;
 
116
      const double kEps = c == 0 ? /* central */ 3e-9 : /* forward */ 2e-5;
 
117
 
 
118
      ExpectClose(x2[i],                    dydx1[5 * 0 + i], kEps);  // y1
 
119
      ExpectClose(x1[i],                    dydx2[5 * 0 + i], kEps);
 
120
      ExpectClose(2 * x2[i] * residuals[0], dydx1[5 * 1 + i], kEps);  // y2
 
121
      ExpectClose(2 * x1[i] * residuals[0], dydx2[5 * 1 + i], kEps);
 
122
      ExpectClose(0.0,                      dydx1[5 * 2 + i], kEps);  // y3
 
123
      ExpectClose(2 * x2[i],                dydx2[5 * 2 + i], kEps);
 
124
    }
 
125
  }
 
126
}
 
127
 
 
128
// y1 = sin(x1'x2)
 
129
// y2 = exp(-x1'x2 / 10)
 
130
//
 
131
// dy1/dx1 =  x2 * cos(x1'x2),            dy1/dx2 =  x1 * cos(x1'x2)
 
132
// dy2/dx1 = -x2 * exp(-x1'x2 / 10) / 10, dy2/dx2 = -x2 * exp(-x1'x2 / 10) / 10
 
133
class TranscendentalTestCostFunction : public CostFunction {
 
134
 public:
 
135
  TranscendentalTestCostFunction() {
 
136
    set_num_residuals(2);
 
137
    mutable_parameter_block_sizes()->push_back(5);  // x1.
 
138
    mutable_parameter_block_sizes()->push_back(5);  // x2.
 
139
  }
 
140
  virtual bool Evaluate(double const* const* parameters,
 
141
                        double* residuals,
 
142
                        double** jacobians) const {
 
143
    (void) jacobians;  // Ignored.
 
144
 
 
145
    double x1x2 = 0;
 
146
    for (int i = 0; i < 5; ++i) {
 
147
      x1x2 += parameters[0][i] * parameters[1][i];
 
148
    }
 
149
    residuals[0] = sin(x1x2);
 
150
    residuals[1] = exp(-x1x2 / 10);
 
151
    return true;
 
152
  }
 
153
};
 
154
 
 
155
TEST(NumericDiffCostFunction, TransendentalOperationsInCostFunction) {
 
156
  // Try both central and forward difference.
 
157
  internal::scoped_ptr<CostFunction> cfs[2];
 
158
  cfs[0].reset(
 
159
      new NumericDiffCostFunction<TranscendentalTestCostFunction,
 
160
                                  CENTRAL,
 
161
                                  2,  /* number of residuals */
 
162
                                  5,  /* size of x1 */
 
163
                                  5   /* size of x2 */>(
 
164
          new TranscendentalTestCostFunction, TAKE_OWNERSHIP));
 
165
 
 
166
  cfs[1].reset(
 
167
      new NumericDiffCostFunction<TranscendentalTestCostFunction,
 
168
                                  FORWARD,
 
169
                                  2,  /* number of residuals */
 
170
                                  5,  /* size of x1 */
 
171
                                  5   /* size of x2 */>(
 
172
          new TranscendentalTestCostFunction, TAKE_OWNERSHIP));
 
173
 
 
174
  for (int c = 0; c < 2; ++c) {
 
175
    CostFunction *cost_function = cfs[c].get();
 
176
 
 
177
    struct {
 
178
      double x1[5];
 
179
      double x2[5];
 
180
    } kTests[] = {
 
181
      { { 1.0, 2.0, 3.0, 4.0, 5.0 },  // No zeros.
 
182
        { 9.0, 9.0, 5.0, 5.0, 1.0 },
 
183
      },
 
184
      { { 0.0, 2.0, 3.0, 0.0, 5.0 },  // Some zeros x1.
 
185
        { 9.0, 9.0, 5.0, 5.0, 1.0 },
 
186
      },
 
187
      { { 1.0, 2.0, 3.0, 1.0, 5.0 },  // Some zeros x2.
 
188
        { 0.0, 9.0, 0.0, 5.0, 0.0 },
 
189
      },
 
190
      { { 0.0, 0.0, 0.0, 0.0, 0.0 },  // All zeros x1.
 
191
        { 9.0, 9.0, 5.0, 5.0, 1.0 },
 
192
      },
 
193
      { { 1.0, 2.0, 3.0, 4.0, 5.0 },  // All zeros x2.
 
194
        { 0.0, 0.0, 0.0, 0.0, 0.0 },
 
195
      },
 
196
      { { 0.0, 0.0, 0.0, 0.0, 0.0 },  // All zeros.
 
197
        { 0.0, 0.0, 0.0, 0.0, 0.0 },
 
198
      },
 
199
    };
 
200
    for (int k = 0; k < ARRAYSIZE(kTests); ++k) {
 
201
      double *x1 = &(kTests[k].x1[0]);
 
202
      double *x2 = &(kTests[k].x2[0]);
 
203
      double *parameters[] = { x1, x2 };
 
204
 
 
205
      double dydx1[10];
 
206
      double dydx2[10];
 
207
      double *jacobians[2] = { &dydx1[0], &dydx2[0] };
 
208
 
 
209
      double residuals[2];
 
210
 
 
211
      ASSERT_TRUE(cost_function->Evaluate(&parameters[0],
 
212
                                          &residuals[0],
 
213
                                          &jacobians[0]));
 
214
      LOG(INFO) << "Ran evaluate for test k=" << k << " c=" << c;
 
215
 
 
216
      double x1x2 = 0;
 
217
      for (int i = 0; i < 5; ++i) {
 
218
        x1x2 += x1[i] * x2[i];
 
219
      }
 
220
 
 
221
      for (int i = 0; i < 5; ++i) {
 
222
        const double kEps = c == 0 ? /* central */ 3e-9 : /* forward */ 2e-5;
 
223
 
 
224
        ExpectClose( x2[i] * cos(x1x2),              dydx1[5 * 0 + i], kEps);
 
225
        ExpectClose( x1[i] * cos(x1x2),              dydx2[5 * 0 + i], kEps);
 
226
        ExpectClose(-x2[i] * exp(-x1x2 / 10.) / 10., dydx1[5 * 1 + i], kEps);
 
227
        ExpectClose(-x1[i] * exp(-x1x2 / 10.) / 10., dydx2[5 * 1 + i], kEps);
 
228
      }
 
229
    }
 
230
  }
 
231
}
 
232
 
 
233
}  // namespace internal
 
234
}  // namespace ceres