Lazily initialize the bounds arrays in ParameterBlock.
Problems that do not use bounds do not have to pay the
price of storing bounds constraints.
Also replace the raw pointer access to the upper and
lower bounds arrays with accessors which hides the
lazy initialization from the user.
Change-Id: I0325a35de9c29f853559f891e32e7c777686e537
diff --git a/internal/ceres/parameter_block.h b/internal/ceres/parameter_block.h
index a978eb5..7bc823d 100644
--- a/internal/ceres/parameter_block.h
+++ b/internal/ceres/parameter_block.h
@@ -184,11 +184,27 @@
void SetUpperBound(int index, double upper_bound) {
CHECK_LT(index, size_);
+
+ if (upper_bounds_.get() == NULL) {
+ upper_bounds_.reset(new double[size_]);
+ std::fill(upper_bounds_.get(),
+ upper_bounds_.get() + size_,
+ std::numeric_limits<double>::max());
+ }
+
upper_bounds_[index] = upper_bound;
};
void SetLowerBound(int index, double lower_bound) {
CHECK_LT(index, size_);
+
+ if (lower_bounds_.get() == NULL) {
+ lower_bounds_.reset(new double[size_]);
+ std::fill(lower_bounds_.get(),
+ lower_bounds_.get() + size_,
+ -std::numeric_limits<double>::max());
+ }
+
lower_bounds_[index] = lower_bound;
}
@@ -206,10 +222,16 @@
}
// Project onto the box constraints.
- for (int i = 0; i < size_; ++i) {
- x_plus_delta[i] = std::min(std::max(x_plus_delta[i],
- lower_bounds_[i]),
- upper_bounds_[i]);
+ if (lower_bounds_.get() != NULL) {
+ for (int i = 0; i < size_; ++i) {
+ x_plus_delta[i] = std::max(x_plus_delta[i], lower_bounds_[i]);
+ }
+ }
+
+ if (upper_bounds_.get() != NULL) {
+ for (int i = 0; i < size_; ++i) {
+ x_plus_delta[i] = std::min(x_plus_delta[i], upper_bounds_[i]);
+ }
}
return true;
@@ -257,12 +279,20 @@
return residual_blocks_.get();
}
- const double* upper_bounds() const {
- return upper_bounds_.get();
+ double LowerBoundForParameter(int index) const {
+ if (lower_bounds_.get() == NULL) {
+ return -std::numeric_limits<double>::max();
+ } else {
+ return lower_bounds_[index];
+ }
}
- const double* lower_bounds() const {
- return lower_bounds_.get();
+ double UpperBoundForParameter(int index) const {
+ if (upper_bounds_.get() == NULL) {
+ return std::numeric_limits<double>::max();
+ } else {
+ return upper_bounds_[index];
+ }
}
private:
@@ -281,15 +311,6 @@
SetParameterization(local_parameterization);
}
- upper_bounds_.reset(new double[size_]);
- std::fill(upper_bounds_.get(),
- upper_bounds_.get() + size_,
- std::numeric_limits<double>::max());
- lower_bounds_.reset(new double[size_]);
- std::fill(lower_bounds_.get(),
- lower_bounds_.get() + size_,
- -std::numeric_limits<double>::max());
-
state_offset_ = -1;
delta_offset_ = -1;
}
@@ -352,8 +373,15 @@
// If non-null, contains the residual blocks this parameter block is in.
scoped_ptr<ResidualBlockSet> residual_blocks_;
- // Upper and lower bounds for the parameter block. These arrays are
- // initialized to std::numeric_limits<double>::max() and
+ // Upper and lower bounds for the parameter block. SetUpperBound
+ // and SetLowerBound lazily initialize the upper_bounds_ and
+ // lower_bounds_ arrays. If they are never called, then memory for
+ // these arrays is never allocated. Thus for problems where there
+ // are no bounds, or only one sided bounds we do not pay the cost of
+ // allocating memory for the inactive bounds constraints.
+ //
+ // Upon initialization these arrays are initialized to
+ // std::numeric_limits<double>::max() and
// -std::numeric_limits<double>::max() respectively which correspond
// to the parameter block being unconstrained.
scoped_array<double> upper_bounds_;
diff --git a/internal/ceres/parameter_block_test.cc b/internal/ceres/parameter_block_test.cc
index a76c408..5a2db3c 100644
--- a/internal/ceres/parameter_block_test.cc
+++ b/internal/ceres/parameter_block_test.cc
@@ -172,26 +172,29 @@
TEST(ParameterBlock, DefaultBounds) {
double x[2];
ParameterBlock parameter_block(x, 2, -1, NULL);
- const double* upper_bounds = parameter_block.upper_bounds();
- EXPECT_EQ(upper_bounds[0], std::numeric_limits<double>::max());
- EXPECT_EQ(upper_bounds[1], std::numeric_limits<double>::max());
- const double* lower_bounds = parameter_block.lower_bounds();
- EXPECT_EQ(lower_bounds[0], -std::numeric_limits<double>::max());
- EXPECT_EQ(lower_bounds[1], -std::numeric_limits<double>::max());
+ EXPECT_EQ(parameter_block.UpperBoundForParameter(0),
+ std::numeric_limits<double>::max());
+ EXPECT_EQ(parameter_block.UpperBoundForParameter(1),
+ std::numeric_limits<double>::max());
+ EXPECT_EQ(parameter_block.LowerBoundForParameter(0),
+ -std::numeric_limits<double>::max());
+ EXPECT_EQ(parameter_block.LowerBoundForParameter(1),
+ -std::numeric_limits<double>::max());
}
TEST(ParameterBlock, SetBounds) {
double x[2];
ParameterBlock parameter_block(x, 2, -1, NULL);
- parameter_block.SetUpperBound(1, 1);
parameter_block.SetLowerBound(0, 1);
+ parameter_block.SetUpperBound(1, 1);
- const double* upper_bounds = parameter_block.upper_bounds();
- EXPECT_EQ(upper_bounds[0], std::numeric_limits<double>::max());
- EXPECT_EQ(upper_bounds[1], 1.0);
- const double* lower_bounds = parameter_block.lower_bounds();
- EXPECT_EQ(lower_bounds[0], 1.0);
- EXPECT_EQ(lower_bounds[1], -std::numeric_limits<double>::max());
+ EXPECT_EQ(parameter_block.LowerBoundForParameter(0), 1.0);
+ EXPECT_EQ(parameter_block.LowerBoundForParameter(1),
+ -std::numeric_limits<double>::max());
+
+ EXPECT_EQ(parameter_block.UpperBoundForParameter(0),
+ std::numeric_limits<double>::max());
+ EXPECT_EQ(parameter_block.UpperBoundForParameter(1), 1.0);
}
TEST(ParameterBlock, PlusWithBoundsConstraints) {
diff --git a/internal/ceres/solver_impl.cc b/internal/ceres/solver_impl.cc
index 6ae2c90..2bf6cd2 100644
--- a/internal/ceres/solver_impl.cc
+++ b/internal/ceres/solver_impl.cc
@@ -328,16 +328,16 @@
bool IsBoundsConstrained(const Program& program) {
const vector<ParameterBlock*>& parameter_blocks = program.parameter_blocks();
for (int i = 0; i < parameter_blocks.size(); ++i) {
- if (parameter_blocks[i]->IsConstant()) {
+ const ParameterBlock* parameter_block = parameter_blocks[i];
+ if (parameter_block->IsConstant()) {
continue;
}
-
- const double* lower_bounds = parameter_blocks[i]->lower_bounds();
- const double* upper_bounds = parameter_blocks[i]->upper_bounds();
- const int size = parameter_blocks[i]->Size();
+ const int size = parameter_block->Size();
for (int j = 0; j < size; ++j) {
- if (lower_bounds[j] > -std::numeric_limits<double>::max() ||
- upper_bounds[j] < std::numeric_limits<double>::max()) {
+ const double lower_bound = parameter_block->LowerBoundForParameter(j);
+ const double upper_bound = parameter_block->UpperBoundForParameter(j);
+ if (lower_bound > -std::numeric_limits<double>::max() ||
+ upper_bound < std::numeric_limits<double>::max()) {
return true;
}
}
@@ -353,24 +353,25 @@
const Program& program = problem->program();
const vector<ParameterBlock*>& parameter_blocks = program.parameter_blocks();
for (int i = 0; i < parameter_blocks.size(); ++i) {
- const double* array = parameter_blocks[i]->user_state();
- const double* lower_bounds = parameter_blocks[i]->lower_bounds();
- const double* upper_bounds = parameter_blocks[i]->upper_bounds();
- const int size = parameter_blocks[i]->Size();
- if (parameter_blocks[i]->IsConstant()) {
+ const ParameterBlock* parameter_block = parameter_blocks[i];
+ const double* parameters = parameter_block->user_state();
+ const int size = parameter_block->Size();
+ if (parameter_block->IsConstant()) {
// Constant parameter blocks must start in the feasible region
// to ultimately produce a feasible solution, since Ceres cannot
// change them.
for (int j = 0; j < size; ++j) {
- if (array[j] < lower_bounds[j] || array[j] > upper_bounds[j]) {
+ const double lower_bound = parameter_block->LowerBoundForParameter(j);
+ const double upper_bound = parameter_block->UpperBoundForParameter(j);
+ if (parameters[j] < lower_bound || parameters[j] > upper_bound) {
*message = StringPrintf(
"ParameterBlock: %p with size %d has at least one infeasible "
"value."
"\nFirst infeasible value is at index: %d."
"\nLower bound: %e, value: %e, upper bound: %e"
"\nParameter block values: ",
- array, size, j, lower_bounds[j], array[j], upper_bounds[j]);
- AppendArrayToString(size, array, message);
+ parameters, size, j, lower_bound, parameters[j], upper_bound);
+ AppendArrayToString(size, parameters, message);
return false;
}
}
@@ -379,15 +380,17 @@
// regions, otherwise there is no way to produce a feasible
// solution.
for (int j = 0; j < size; ++j) {
- if (lower_bounds[j] >= upper_bounds[j]) {
+ const double lower_bound = parameter_block->LowerBoundForParameter(j);
+ const double upper_bound = parameter_block->UpperBoundForParameter(j);
+ if (lower_bound >= upper_bound) {
*message = StringPrintf(
"ParameterBlock: %p with size %d has at least one infeasible "
"bound."
"\nFirst infeasible bound is at index: %d."
"\nLower bound: %e, upper bound: %e"
"\nParameter block values: ",
- array, size, j, lower_bounds[j], upper_bounds[j]);
- AppendArrayToString(size, array, message);
+ parameters, size, j, lower_bound, upper_bound);
+ AppendArrayToString(size, parameters, message);
return false;
}
}
