Reduce the number of individual PRNG instances
Use same instance of a PRNG throughout by passing it to methods and
functions as an argument to generate random numbers without breaking the
sequence.
Change-Id: Ib024bbc1ea2d14e4b9afb71857856a5fb77b1667
diff --git a/internal/ceres/autodiff_test.cc b/internal/ceres/autodiff_test.cc
index 629eea8..0916a45 100644
--- a/internal/ceres/autodiff_test.cc
+++ b/internal/ceres/autodiff_test.cc
@@ -30,6 +30,8 @@
#include "ceres/internal/autodiff.h"
+#include <algorithm>
+#include <iterator>
#include <random>
#include "gtest/gtest.h"
@@ -173,9 +175,9 @@
// Make random P and X, in a single vector.
double PX[12 + 4];
- for (double& PX_i : PX) {
- PX_i = uniform01(prng);
- }
+ std::generate(std::begin(PX), std::end(PX), [&prng, &uniform01] {
+ return uniform01(prng);
+ });
// Handy names for the P and X parts.
double* P = PX + 0;
@@ -294,7 +296,10 @@
double qcX[4 + 3 + 3];
std::mt19937 prng;
std::uniform_real_distribution<double> uniform01(0.0, 1.0);
- for (double& qcX_i : qcX) qcX_i = uniform01(prng);
+
+ std::generate(std::begin(qcX), std::end(qcX), [&prng, &uniform01] {
+ return uniform01(prng);
+ });
// Handy names.
double* q = qcX;
diff --git a/internal/ceres/block_sparse_matrix.cc b/internal/ceres/block_sparse_matrix.cc
index ae6bd3a..90bd926 100644
--- a/internal/ceres/block_sparse_matrix.cc
+++ b/internal/ceres/block_sparse_matrix.cc
@@ -375,7 +375,7 @@
}
std::unique_ptr<BlockSparseMatrix> BlockSparseMatrix::CreateRandomMatrix(
- const BlockSparseMatrix::RandomMatrixOptions& options) {
+ const BlockSparseMatrix::RandomMatrixOptions& options, std::mt19937& prng) {
CHECK_GT(options.num_row_blocks, 0);
CHECK_GT(options.min_row_block_size, 0);
CHECK_GT(options.max_row_block_size, 0);
@@ -383,7 +383,6 @@
CHECK_GT(options.block_density, 0.0);
CHECK_LE(options.block_density, 1.0);
- std::mt19937 prng;
std::uniform_int_distribution<int> col_distribution(
options.min_col_block_size, options.max_col_block_size);
std::uniform_int_distribution<int> row_distribution(
@@ -436,9 +435,10 @@
auto matrix = std::make_unique<BlockSparseMatrix>(bs);
double* values = matrix->mutable_values();
std::normal_distribution<double> standard_normal_distribution;
- for (int i = 0; i < matrix->num_nonzeros(); ++i) {
- values[i] = standard_normal_distribution(prng);
- }
+ std::generate_n(
+ values, matrix->num_nonzeros(), [&standard_normal_distribution, &prng] {
+ return standard_normal_distribution(prng);
+ });
return matrix;
}
diff --git a/internal/ceres/block_sparse_matrix.h b/internal/ceres/block_sparse_matrix.h
index 7cef18d..44cb5e5 100644
--- a/internal/ceres/block_sparse_matrix.h
+++ b/internal/ceres/block_sparse_matrix.h
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -35,6 +35,7 @@
#define CERES_INTERNAL_BLOCK_SPARSE_MATRIX_H_
#include <memory>
+#include <random>
#include "ceres/block_structure.h"
#include "ceres/crs_matrix.h"
@@ -123,7 +124,7 @@
// distributed and whose structure is determined by
// RandomMatrixOptions.
static std::unique_ptr<BlockSparseMatrix> CreateRandomMatrix(
- const RandomMatrixOptions& options);
+ const RandomMatrixOptions& options, std::mt19937& prng);
private:
int num_rows_;
diff --git a/internal/ceres/compressed_row_sparse_matrix.cc b/internal/ceres/compressed_row_sparse_matrix.cc
index 574b0c0..ce587da 100644
--- a/internal/ceres/compressed_row_sparse_matrix.cc
+++ b/internal/ceres/compressed_row_sparse_matrix.cc
@@ -119,11 +119,12 @@
transpose_rows[0] = 0;
}
+template <class RandomNormalFunctor>
void AddRandomBlock(const int num_rows,
const int num_cols,
const int row_block_begin,
const int col_block_begin,
- std::function<double()> dist,
+ RandomNormalFunctor&& randn,
std::vector<int>* rows,
std::vector<int>* cols,
std::vector<double>* values) {
@@ -131,20 +132,21 @@
for (int c = 0; c < num_cols; ++c) {
rows->push_back(row_block_begin + r);
cols->push_back(col_block_begin + c);
- values->push_back(dist());
+ values->push_back(randn());
}
}
}
+template <class RandomNormalFunctor>
void AddSymmetricRandomBlock(const int num_rows,
const int row_block_begin,
- std::function<double()> dist,
+ RandomNormalFunctor&& randn,
std::vector<int>* rows,
std::vector<int>* cols,
std::vector<double>* values) {
for (int r = 0; r < num_rows; ++r) {
for (int c = r; c < num_rows; ++c) {
- const double v = dist();
+ const double v = randn();
rows->push_back(row_block_begin + r);
cols->push_back(row_block_begin + c);
values->push_back(v);
@@ -625,7 +627,8 @@
std::unique_ptr<CompressedRowSparseMatrix>
CompressedRowSparseMatrix::CreateRandomMatrix(
- CompressedRowSparseMatrix::RandomMatrixOptions options) {
+ CompressedRowSparseMatrix::RandomMatrixOptions options,
+ std::mt19937& prng) {
CHECK_GT(options.num_row_blocks, 0);
CHECK_GT(options.min_row_block_size, 0);
CHECK_GT(options.max_row_block_size, 0);
@@ -651,7 +654,6 @@
vector<int> col_blocks;
col_blocks.reserve(options.num_col_blocks);
- std::mt19937 prng;
std::uniform_int_distribution<int> col_distribution(
options.min_col_block_size, options.max_col_block_size);
std::uniform_int_distribution<int> row_distribution(
@@ -708,6 +710,9 @@
// Randomly determine if this block is present or not.
if (uniform01(prng) <= options.block_density) {
+ auto randn = [&standard_normal_distribution, &prng] {
+ return standard_normal_distribution(prng);
+ };
// If the matrix is symmetric, then we take care to generate
// symmetric diagonal blocks.
if (options.storage_type == StorageType::UNSYMMETRIC || r != c) {
@@ -715,14 +720,14 @@
col_blocks[c],
row_block_begin,
col_block_begin,
- values_dist,
+ randn,
&tsm_rows,
&tsm_cols,
&tsm_values);
} else {
AddSymmetricRandomBlock(row_blocks[r],
row_block_begin,
- values_dist,
+ randn,
&tsm_rows,
&tsm_cols,
&tsm_values);
diff --git a/internal/ceres/compressed_row_sparse_matrix.h b/internal/ceres/compressed_row_sparse_matrix.h
index 2580045..d164bef 100644
--- a/internal/ceres/compressed_row_sparse_matrix.h
+++ b/internal/ceres/compressed_row_sparse_matrix.h
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -32,6 +32,7 @@
#define CERES_INTERNAL_COMPRESSED_ROW_SPARSE_MATRIX_H_
#include <memory>
+#include <random>
#include <vector>
#include "ceres/internal/disable_warnings.h"
@@ -195,7 +196,7 @@
// normally distributed and whose structure is determined by
// RandomMatrixOptions.
static std::unique_ptr<CompressedRowSparseMatrix> CreateRandomMatrix(
- RandomMatrixOptions options);
+ RandomMatrixOptions options, std::mt19937& prng);
private:
static std::unique_ptr<CompressedRowSparseMatrix> FromTripletSparseMatrix(
diff --git a/internal/ceres/compressed_row_sparse_matrix_test.cc b/internal/ceres/compressed_row_sparse_matrix_test.cc
index 42f5498..718a659 100644
--- a/internal/ceres/compressed_row_sparse_matrix_test.cc
+++ b/internal/ceres/compressed_row_sparse_matrix_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -328,6 +328,7 @@
}
TEST(CompressedRowSparseMatrix, FromTripletSparseMatrix) {
+ std::mt19937 prng;
TripletSparseMatrix::RandomMatrixOptions options;
options.num_rows = 5;
options.num_cols = 7;
@@ -336,7 +337,7 @@
const int kNumTrials = 10;
for (int i = 0; i < kNumTrials; ++i) {
std::unique_ptr<TripletSparseMatrix> tsm =
- TripletSparseMatrix::CreateRandomMatrix(options);
+ TripletSparseMatrix::CreateRandomMatrix(options, prng);
std::unique_ptr<CompressedRowSparseMatrix> crsm =
CompressedRowSparseMatrix::FromTripletSparseMatrix(*tsm);
@@ -354,6 +355,7 @@
}
TEST(CompressedRowSparseMatrix, FromTripletSparseMatrixTransposed) {
+ std::mt19937 prng;
TripletSparseMatrix::RandomMatrixOptions options;
options.num_rows = 5;
options.num_cols = 7;
@@ -362,7 +364,7 @@
const int kNumTrials = 10;
for (int i = 0; i < kNumTrials; ++i) {
std::unique_ptr<TripletSparseMatrix> tsm =
- TripletSparseMatrix::CreateRandomMatrix(options);
+ TripletSparseMatrix::CreateRandomMatrix(options, prng);
std::unique_ptr<CompressedRowSparseMatrix> crsm =
CompressedRowSparseMatrix::FromTripletSparseMatrixTransposed(*tsm);
@@ -421,7 +423,7 @@
options.block_density = uniform(prng);
options.storage_type = ::testing::get<0>(param);
std::unique_ptr<CompressedRowSparseMatrix> matrix =
- CompressedRowSparseMatrix::CreateRandomMatrix(options);
+ CompressedRowSparseMatrix::CreateRandomMatrix(options, prng);
const int num_rows = matrix->num_rows();
const int num_cols = matrix->num_cols();
@@ -491,7 +493,7 @@
options.block_density = uniform(prng);
options.storage_type = ::testing::get<0>(param);
std::unique_ptr<CompressedRowSparseMatrix> matrix =
- CompressedRowSparseMatrix::CreateRandomMatrix(options);
+ CompressedRowSparseMatrix::CreateRandomMatrix(options, prng);
const int num_rows = matrix->num_rows();
const int num_cols = matrix->num_cols();
@@ -561,7 +563,7 @@
options.block_density = uniform(prng);
options.storage_type = ::testing::get<0>(param);
std::unique_ptr<CompressedRowSparseMatrix> matrix =
- CompressedRowSparseMatrix::CreateRandomMatrix(options);
+ CompressedRowSparseMatrix::CreateRandomMatrix(options, prng);
const int num_cols = matrix->num_cols();
Vector actual(num_cols);
diff --git a/internal/ceres/corrector_test.cc b/internal/ceres/corrector_test.cc
index 7f58efa..a5dd29c 100644
--- a/internal/ceres/corrector_test.cc
+++ b/internal/ceres/corrector_test.cc
@@ -32,7 +32,6 @@
#include <algorithm>
#include <cmath>
-#include <cstdlib>
#include <cstring>
#include <random>
diff --git a/internal/ceres/gradient_checker_test.cc b/internal/ceres/gradient_checker_test.cc
index 3c5f5da..8d7d10f 100644
--- a/internal/ceres/gradient_checker_test.cc
+++ b/internal/ceres/gradient_checker_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2021 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -33,7 +33,6 @@
#include "ceres/gradient_checker.h"
#include <cmath>
-#include <cstdlib>
#include <random>
#include <utility>
#include <vector>
@@ -59,16 +58,16 @@
// version, they are both block vectors, of course.
class GoodTestTerm : public CostFunction {
public:
- GoodTestTerm(int arity, int const* dim) : arity_(arity), return_value_(true) {
- std::mt19937 prng;
- std::uniform_real_distribution<double> distribution(-1.0, 1.0);
-
+ template <class UniformRandomFunctor>
+ GoodTestTerm(int arity, int const* dim, UniformRandomFunctor&& randu)
+ : arity_(arity), return_value_(true) {
+ std::uniform_real_distribution distribution(-1.0, 1.0);
// Make 'arity' random vectors.
a_.resize(arity_);
for (int j = 0; j < arity_; ++j) {
a_[j].resize(dim[j]);
for (int u = 0; u < dim[j]; ++u) {
- a_[j][u] = distribution(prng);
+ a_[j][u] = randu();
}
}
@@ -121,15 +120,15 @@
class BadTestTerm : public CostFunction {
public:
- BadTestTerm(int arity, int const* dim) : arity_(arity) {
- std::mt19937 prng;
- std::uniform_real_distribution<double> distribution(-1.0, 1.0);
+ template <class UniformRandomFunctor>
+ BadTestTerm(int arity, int const* dim, UniformRandomFunctor&& randu)
+ : arity_(arity) {
// Make 'arity' random vectors.
a_.resize(arity_);
for (int j = 0; j < arity_; ++j) {
a_[j].resize(dim[j]);
for (int u = 0; u < dim[j]; ++u) {
- a_[j][u] = distribution(prng);
+ a_[j][u] = randu();
}
}
@@ -210,10 +209,11 @@
FixedArray<double*> parameters(num_parameters);
std::mt19937 prng;
std::uniform_real_distribution<double> distribution(-1.0, 1.0);
+ auto randu = [&prng, &distribution] { return distribution(prng); };
for (int j = 0; j < num_parameters; ++j) {
parameters[j] = new double[parameter_sizes[j]];
for (int u = 0; u < parameter_sizes[j]; ++u) {
- parameters[j][u] = distribution(prng);
+ parameters[j][u] = randu();
}
}
@@ -221,7 +221,7 @@
GradientChecker::ProbeResults results;
// Test that Probe returns true for correct Jacobians.
- GoodTestTerm good_term(num_parameters, parameter_sizes.data());
+ GoodTestTerm good_term(num_parameters, parameter_sizes.data(), randu);
std::vector<const Manifold*>* manifolds = nullptr;
GradientChecker good_gradient_checker(
&good_term, manifolds, numeric_diff_options);
@@ -258,7 +258,7 @@
EXPECT_FALSE(results.error_log.empty());
// Test that Probe returns false for incorrect Jacobians.
- BadTestTerm bad_term(num_parameters, parameter_sizes.data());
+ BadTestTerm bad_term(num_parameters, parameter_sizes.data(), randu);
GradientChecker bad_gradient_checker(
&bad_term, manifolds, numeric_diff_options);
EXPECT_FALSE(
diff --git a/internal/ceres/gradient_checking_cost_function_test.cc b/internal/ceres/gradient_checking_cost_function_test.cc
index caf4f25..87c3278 100644
--- a/internal/ceres/gradient_checking_cost_function_test.cc
+++ b/internal/ceres/gradient_checking_cost_function_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -69,15 +69,15 @@
public:
// The constructor of this function needs to know the number
// of blocks desired, and the size of each block.
- TestTerm(int arity, int const* dim) : arity_(arity) {
- std::mt19937 prng;
- std::uniform_real_distribution<double> uniform01(-1.0, 1.0);
+ template <class UniformRandomFunctor>
+ TestTerm(int arity, int const* dim, UniformRandomFunctor&& randu)
+ : arity_(arity) {
// Make 'arity' random vectors.
a_.resize(arity_);
for (int j = 0; j < arity_; ++j) {
a_[j].resize(dim[j]);
for (int u = 0; u < dim[j]; ++u) {
- a_[j][u] = uniform01(prng);
+ a_[j][u] = randu();
}
}
@@ -139,11 +139,12 @@
// Make a random set of blocks.
vector<double*> parameters(arity);
std::mt19937 prng;
- std::uniform_real_distribution<double> uniform01(-1.0, 1.0);
+ std::uniform_real_distribution<double> distribution(-1.0, 1.0);
+ auto randu = [&prng, &distribution] { return distribution(prng); };
for (int j = 0; j < arity; ++j) {
parameters[j] = new double[dim[j]];
for (int u = 0; u < dim[j]; ++u) {
- parameters[j][u] = uniform01(prng);
+ parameters[j][u] = randu();
}
}
@@ -162,7 +163,7 @@
const double kRelativeStepSize = 1e-6;
const double kRelativePrecision = 1e-4;
- TestTerm<-1, -1> term(arity, dim);
+ TestTerm<-1, -1> term(arity, dim, randu);
GradientCheckingIterationCallback callback;
auto gradient_checking_cost_function =
CreateGradientCheckingCostFunction(&term,
@@ -196,11 +197,12 @@
// Make a random set of blocks.
vector<double*> parameters(arity);
std::mt19937 prng;
- std::uniform_real_distribution<double> uniform01(-1.0, 1.0);
+ std::uniform_real_distribution<double> distribution(-1.0, 1.0);
+ auto randu = [&prng, &distribution] { return distribution(prng); };
for (int j = 0; j < arity; ++j) {
parameters[j] = new double[dim[j]];
for (int u = 0; u < dim[j]; ++u) {
- parameters[j][u] = uniform01(prng);
+ parameters[j][u] = randu();
}
}
@@ -218,7 +220,7 @@
// Should have one term that's bad, causing everything to get dumped.
LOG(INFO) << "Bad gradient";
{
- TestTerm<1, 2> term(arity, dim);
+ TestTerm<1, 2> term(arity, dim, randu);
GradientCheckingIterationCallback callback;
auto gradient_checking_cost_function =
CreateGradientCheckingCostFunction(&term,
@@ -238,7 +240,7 @@
// The gradient is correct, so no errors are reported.
LOG(INFO) << "Good gradient";
{
- TestTerm<-1, -1> term(arity, dim);
+ TestTerm<-1, -1> term(arity, dim, randu);
GradientCheckingIterationCallback callback;
auto gradient_checking_cost_function =
CreateGradientCheckingCostFunction(&term,
diff --git a/internal/ceres/inner_product_computer_test.cc b/internal/ceres/inner_product_computer_test.cc
index 5bd4c50..4573b6b 100644
--- a/internal/ceres/inner_product_computer_test.cc
+++ b/internal/ceres/inner_product_computer_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2017 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -109,7 +109,7 @@
VLOG(2) << "block density: " << options.block_density;
std::unique_ptr<BlockSparseMatrix> random_matrix(
- BlockSparseMatrix::CreateRandomMatrix(options));
+ BlockSparseMatrix::CreateRandomMatrix(options, prng));
TripletSparseMatrix tsm(random_matrix->num_rows(),
random_matrix->num_cols(),
@@ -169,7 +169,7 @@
VLOG(2) << "block density: " << options.block_density;
std::unique_ptr<BlockSparseMatrix> random_matrix(
- BlockSparseMatrix::CreateRandomMatrix(options));
+ BlockSparseMatrix::CreateRandomMatrix(options, prng));
const std::vector<CompressedRow>& row_blocks =
random_matrix->block_structure()->rows;
diff --git a/internal/ceres/normal_prior_test.cc b/internal/ceres/normal_prior_test.cc
index 25e6322..881038f 100644
--- a/internal/ceres/normal_prior_test.cc
+++ b/internal/ceres/normal_prior_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -42,6 +42,7 @@
TEST(NormalPriorTest, ResidualAtRandomPosition) {
std::mt19937 prng;
std::uniform_real_distribution<double> distribution(-1.0, 1.0);
+ auto randu = [&distribution, &prng] { return distribution(prng); };
for (int num_rows = 1; num_rows < 5; ++num_rows) {
for (int num_cols = 1; num_cols < 5; ++num_cols) {
Vector b(num_cols);
@@ -50,7 +51,7 @@
A.setRandom();
auto* x = new double[num_cols];
- for (int i = 0; i < num_cols; ++i) x[i] = distribution(prng);
+ std::generate_n(x, num_cols, randu);
auto* jacobian = new double[num_rows * num_cols];
Vector residuals(num_rows);
@@ -77,6 +78,7 @@
TEST(NormalPriorTest, ResidualAtRandomPositionNullJacobians) {
std::mt19937 prng;
std::uniform_real_distribution<double> distribution(-1.0, 1.0);
+ auto randu = [&distribution, &prng] { return distribution(prng); };
for (int num_rows = 1; num_rows < 5; ++num_rows) {
for (int num_cols = 1; num_cols < 5; ++num_cols) {
Vector b(num_cols);
@@ -85,7 +87,7 @@
A.setRandom();
auto* x = new double[num_cols];
- for (int i = 0; i < num_cols; ++i) x[i] = distribution(prng);
+ std::generate_n(x, num_cols, randu);
double* jacobians[1];
jacobians[0] = nullptr;
diff --git a/internal/ceres/numeric_diff_cost_function_test.cc b/internal/ceres/numeric_diff_cost_function_test.cc
index b9a0beb..dcffcff 100644
--- a/internal/ceres/numeric_diff_cost_function_test.cc
+++ b/internal/ceres/numeric_diff_cost_function_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -35,6 +35,7 @@
#include <array>
#include <cmath>
#include <memory>
+#include <random>
#include <string>
#include <vector>
@@ -310,6 +311,7 @@
}
TEST(NumericDiffCostFunction, RandomizedFunctorRidders) {
+ std::mt19937 prng;
NumericDiffOptions options;
// Larger initial step size is chosen to produce robust results in the
// presence of random noise.
@@ -321,15 +323,16 @@
1, // number of residuals
1 // size of x1
>>(
- new RandomizedFunctor(kNoiseFactor, kRandomSeed),
+ new RandomizedFunctor(kNoiseFactor, prng),
TAKE_OWNERSHIP,
1,
options);
- RandomizedFunctor functor(kNoiseFactor, kRandomSeed);
+ RandomizedFunctor functor(kNoiseFactor, prng);
functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function);
}
TEST(NumericDiffCostFunction, RandomizedCostFunctionRidders) {
+ std::mt19937 prng;
NumericDiffOptions options;
// Larger initial step size is chosen to produce robust results in the
// presence of random noise.
@@ -341,12 +344,12 @@
1, // number of residuals
1 // size of x1
>>(
- new RandomizedCostFunction(kNoiseFactor, kRandomSeed),
+ new RandomizedCostFunction(kNoiseFactor, prng),
TAKE_OWNERSHIP,
1,
options);
- RandomizedFunctor functor(kNoiseFactor, kRandomSeed);
+ RandomizedFunctor functor(kNoiseFactor, prng);
functor.ExpectCostFunctionEvaluationIsNearlyCorrect(*cost_function);
}
diff --git a/internal/ceres/numeric_diff_test_utils.cc b/internal/ceres/numeric_diff_test_utils.cc
index 3d78e36..9da514e 100644
--- a/internal/ceres/numeric_diff_test_utils.cc
+++ b/internal/ceres/numeric_diff_test_utils.cc
@@ -226,14 +226,7 @@
}
bool RandomizedFunctor::operator()(const double* x1, double* residuals) const {
- double random_value =
- static_cast<double>(rand()) / static_cast<double>(RAND_MAX);
-
- // Normalize noise to [-factor, factor].
- random_value *= 2.0;
- random_value -= 1.0;
- random_value *= noise_factor_;
-
+ double random_value = uniform_distribution_(*prng_);
residuals[0] = x1[0] * x1[0] + random_value;
return true;
}
@@ -244,9 +237,6 @@
const double kTolerance = 2e-4;
- // Initialize random number generator with given seed.
- srand(random_seed_);
-
for (double& test : kTests) {
double* parameters[] = {&test};
double dydx;
diff --git a/internal/ceres/numeric_diff_test_utils.h b/internal/ceres/numeric_diff_test_utils.h
index 68ae5ed..8d3f14b 100644
--- a/internal/ceres/numeric_diff_test_utils.h
+++ b/internal/ceres/numeric_diff_test_utils.h
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -31,6 +31,8 @@
#ifndef CERES_INTERNAL_NUMERIC_DIFF_TEST_UTILS_H_
#define CERES_INTERNAL_NUMERIC_DIFF_TEST_UTILS_H_
+#include <random>
+
#include "ceres/cost_function.h"
#include "ceres/internal/export.h"
#include "ceres/sized_cost_function.h"
@@ -116,8 +118,10 @@
// y = x^2 + [random noise], dy/dx ~ 2x
class CERES_NO_EXPORT RandomizedFunctor {
public:
- RandomizedFunctor(double noise_factor, unsigned int random_seed)
- : noise_factor_(noise_factor), random_seed_(random_seed) {}
+ RandomizedFunctor(double noise_factor, std::mt19937& prng)
+ : noise_factor_(noise_factor),
+ prng_(&prng),
+ uniform_distribution_{-noise_factor, noise_factor} {}
bool operator()(const double* x1, double* residuals) const;
void ExpectCostFunctionEvaluationIsNearlyCorrect(
@@ -125,13 +129,16 @@
private:
double noise_factor_;
- unsigned int random_seed_;
+ // Store the generator as a pointer to be able to modify the instance the
+ // pointer is pointing to.
+ std::mt19937* prng_;
+ mutable std::uniform_real_distribution<> uniform_distribution_;
};
class CERES_NO_EXPORT RandomizedCostFunction : public SizedCostFunction<1, 1> {
public:
- RandomizedCostFunction(double noise_factor, unsigned int random_seed)
- : functor_(noise_factor, random_seed) {}
+ RandomizedCostFunction(double noise_factor, std::mt19937& prng)
+ : functor_(noise_factor, prng) {}
bool Evaluate(double const* const* parameters,
double* residuals,
diff --git a/internal/ceres/partitioned_matrix_view_test.cc b/internal/ceres/partitioned_matrix_view_test.cc
index cb3dd14..e4484dd 100644
--- a/internal/ceres/partitioned_matrix_view_test.cc
+++ b/internal/ceres/partitioned_matrix_view_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
diff --git a/internal/ceres/reorder_program_test.cc b/internal/ceres/reorder_program_test.cc
index 872f381..d6a6577 100644
--- a/internal/ceres/reorder_program_test.cc
+++ b/internal/ceres/reorder_program_test.cc
@@ -265,7 +265,7 @@
problem.GetResidualBlocks(&residual_block_ids);
std::vector<ResidualBlock*> residual_blocks =
problem.program().residual_blocks();
- auto rng = std::default_random_engine{};
+ auto rng = std::mt19937{};
for (int i = 1; i < 6; ++i) {
std::shuffle(
std::begin(residual_block_ids), std::end(residual_block_ids), rng);
diff --git a/internal/ceres/rotation_test.cc b/internal/ceres/rotation_test.cc
index d6de650..97ea514 100644
--- a/internal/ceres/rotation_test.cc
+++ b/internal/ceres/rotation_test.cc
@@ -32,6 +32,7 @@
#include <cmath>
#include <limits>
+#include <random>
#include <string>
#include "ceres/internal/eigen.h"
@@ -56,11 +57,6 @@
const double kPi = 3.14159265358979323846;
const double kHalfSqrt2 = 0.707106781186547524401;
-static double RandDouble() {
- double r = rand();
- return r / RAND_MAX;
-}
-
// A tolerance value for floating-point comparisons.
static double const kTolerance = numeric_limits<double>::epsilon() * 10;
@@ -322,20 +318,22 @@
// Takes a bunch of random axis/angle values, converts them to quaternions,
// and back again.
TEST(Rotation, AngleAxisToQuaterionAndBack) {
- srand(5);
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
for (int i = 0; i < kNumTrials; i++) {
double axis_angle[3];
// Make an axis by choosing three random numbers in [-1, 1) and
// normalizing.
double norm = 0;
for (double& coeff : axis_angle) {
- coeff = RandDouble() * 2 - 1;
+ coeff = uniform_distribution(prng);
norm += coeff * coeff;
}
norm = sqrt(norm);
// Angle in [-pi, pi).
- double theta = kPi * 2 * RandDouble() - kPi;
+ double theta = uniform_distribution(
+ prng, std::uniform_real_distribution<double>::param_type{-kPi, kPi});
for (double& coeff : axis_angle) {
coeff = coeff * theta / norm;
}
@@ -355,13 +353,14 @@
// Takes a bunch of random quaternions, converts them to axis/angle,
// and back again.
TEST(Rotation, QuaterionToAngleAxisAndBack) {
- srand(5);
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
for (int i = 0; i < kNumTrials; i++) {
double quaternion[4];
// Choose four random numbers in [-1, 1) and normalize.
double norm = 0;
for (double& coeff : quaternion) {
- coeff = RandDouble() * 2 - 1;
+ coeff = uniform_distribution(prng);
norm += coeff * coeff;
}
norm = sqrt(norm);
@@ -432,20 +431,24 @@
double matrix[9];
double out_axis_angle[3];
- srand(5);
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
for (int i = 0; i < kNumTrials; i++) {
// Make an axis by choosing three random numbers in [-1, 1) and
// normalizing.
double norm = 0;
for (double& coeff : in_axis_angle) {
- coeff = RandDouble() * 2 - 1;
+ coeff = uniform_distribution(prng);
norm += coeff * coeff;
}
norm = sqrt(norm);
// Angle in [pi - kMaxSmallAngle, pi).
- const double kMaxSmallAngle = 1e-8;
- double theta = kPi - kMaxSmallAngle * RandDouble();
+ constexpr double kMaxSmallAngle = 1e-8;
+ double theta =
+ uniform_distribution(prng,
+ std::uniform_real_distribution<double>::param_type{
+ kPi - kMaxSmallAngle, kPi});
for (double& coeff : in_axis_angle) {
coeff *= (theta / norm);
@@ -531,20 +534,22 @@
// Takes a bunch of random axis/angle values, converts them to rotation
// matrices, and back again.
TEST(Rotation, AngleAxisToRotationMatrixAndBack) {
- srand(5);
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
for (int i = 0; i < kNumTrials; i++) {
double axis_angle[3];
// Make an axis by choosing three random numbers in [-1, 1) and
// normalizing.
double norm = 0;
for (double& i : axis_angle) {
- i = RandDouble() * 2 - 1;
+ i = uniform_distribution(prng);
norm += i * i;
}
norm = sqrt(norm);
// Angle in [-pi, pi).
- double theta = kPi * 2 * RandDouble() - kPi;
+ double theta = uniform_distribution(
+ prng, std::uniform_real_distribution<double>::param_type{-kPi, kPi});
for (double& i : axis_angle) {
i = i * theta / norm;
}
@@ -564,20 +569,25 @@
// Takes a bunch of random axis/angle values near zero, converts them
// to rotation matrices, and back again.
TEST(Rotation, AngleAxisToRotationMatrixAndBackNearZero) {
- srand(5);
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
for (int i = 0; i < kNumTrials; i++) {
double axis_angle[3];
// Make an axis by choosing three random numbers in [-1, 1) and
// normalizing.
double norm = 0;
for (double& i : axis_angle) {
- i = RandDouble() * 2 - 1;
+ i = uniform_distribution(prng);
norm += i * i;
}
norm = sqrt(norm);
// Tiny theta.
- double theta = 1e-16 * (kPi * 2 * RandDouble() - kPi);
+ constexpr double kScale = 1e-16;
+ double theta =
+ uniform_distribution(prng,
+ std::uniform_real_distribution<double>::param_type{
+ -kScale * kPi, kScale * kPi});
for (double& i : axis_angle) {
i = i * theta / norm;
}
@@ -647,11 +657,12 @@
// Test that a random rotation produces an orthonormal rotation
// matrix.
TEST(EulerAnglesToRotationMatrix, IsOrthonormal) {
- srand(5);
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-180.0, 180.0};
for (int trial = 0; trial < kNumTrials; ++trial) {
double euler_angles_degrees[3];
for (double& euler_angles_degree : euler_angles_degrees) {
- euler_angles_degree = RandDouble() * 360.0 - 180.0;
+ euler_angles_degree = uniform_distribution(prng);
}
double rotation_matrix[9];
EulerAnglesToRotationMatrix(euler_angles_degrees, 3, rotation_matrix);
@@ -906,13 +917,15 @@
// Verify that (a * b) * c == a * (b * c).
TEST(Quaternion, MultiplicationIsAssociative) {
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
double a[4];
double b[4];
double c[4];
for (int i = 0; i < 4; ++i) {
- a[i] = 2 * RandDouble() - 1;
- b[i] = 2 * RandDouble() - 1;
- c[i] = 2 * RandDouble() - 1;
+ a[i] = uniform_distribution(prng);
+ b[i] = uniform_distribution(prng);
+ c[i] = uniform_distribution(prng);
}
double ab[4];
@@ -932,6 +945,8 @@
}
TEST(AngleAxis, RotatePointGivesSameAnswerAsRotationMatrix) {
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
double angle_axis[3];
double R[9];
double p[3];
@@ -943,8 +958,8 @@
for (int j = 0; j < 50; ++j) {
double norm2 = 0.0;
for (int k = 0; k < 3; ++k) {
- angle_axis[k] = 2.0 * RandDouble() - 1.0;
- p[k] = 2.0 * RandDouble() - 1.0;
+ angle_axis[k] = uniform_distribution(prng);
+ p[k] = uniform_distribution(prng);
norm2 = angle_axis[k] * angle_axis[k];
}
@@ -976,6 +991,8 @@
}
TEST(AngleAxis, NearZeroRotatePointGivesSameAnswerAsRotationMatrix) {
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution{-1.0, 1.0};
double angle_axis[3];
double R[9];
double p[3];
@@ -985,8 +1002,8 @@
for (int i = 0; i < 10000; ++i) {
double norm2 = 0.0;
for (int k = 0; k < 3; ++k) {
- angle_axis[k] = 2.0 * RandDouble() - 1.0;
- p[k] = 2.0 * RandDouble() - 1.0;
+ angle_axis[k] = uniform_distribution(prng);
+ p[k] = uniform_distribution(prng);
norm2 = angle_axis[k] * angle_axis[k];
}
@@ -1111,7 +1128,12 @@
}
TEST(RotationMatrixToAngleAxis, ExhaustiveRoundTrip) {
- const double kMaxSmallAngle = 1e-8;
+ constexpr double kMaxSmallAngle = 1e-8;
+ std::mt19937 prng;
+ std::uniform_real_distribution<double> uniform_distribution1{
+ kPi - kMaxSmallAngle, kPi};
+ std::uniform_real_distribution<double> uniform_distribution2{
+ -1.0, 2.0 * kMaxSmallAngle - 1.0};
const int kNumSteps = 1000;
for (int i = 0; i < kNumSteps; ++i) {
const double theta = static_cast<double>(i) / kNumSteps * 2.0 * kPi;
@@ -1121,10 +1143,10 @@
CheckRotationMatrixToAngleAxisRoundTrip(theta, phi, kPi);
// Rotation of angle approximately Pi.
CheckRotationMatrixToAngleAxisRoundTrip(
- theta, phi, kPi - kMaxSmallAngle * RandDouble());
+ theta, phi, uniform_distribution1(prng));
// Rotations of angle approximately zero.
CheckRotationMatrixToAngleAxisRoundTrip(
- theta, phi, kMaxSmallAngle * 2.0 * RandDouble() - 1.0);
+ theta, phi, uniform_distribution2(prng));
}
}
}
diff --git a/internal/ceres/schur_eliminator_benchmark.cc b/internal/ceres/schur_eliminator_benchmark.cc
index 0973e73..c6425ce 100644
--- a/internal/ceres/schur_eliminator_benchmark.cc
+++ b/internal/ceres/schur_eliminator_benchmark.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2019 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -91,11 +91,9 @@
matrix_ = std::make_unique<BlockSparseMatrix>(bs);
double* values = matrix_->mutable_values();
- std::mt19937 prng;
- std::normal_distribution<double> standard_normal;
- for (int i = 0; i < matrix_->num_nonzeros(); ++i) {
- values[i] = standard_normal(prng);
- }
+ std::generate_n(values, matrix_->num_nonzeros(), [this] {
+ return standard_normal(prng);
+ });
b_.resize(matrix_->num_rows());
b_.setRandom();
@@ -128,6 +126,8 @@
Vector diagonal_;
Vector z_;
Vector y_;
+ std::mt19937 prng;
+ std::normal_distribution<> standard_normal;
};
static void BM_SchurEliminatorEliminate(benchmark::State& state) {
diff --git a/internal/ceres/schur_eliminator_test.cc b/internal/ceres/schur_eliminator_test.cc
index da27019..17ce67d 100644
--- a/internal/ceres/schur_eliminator_test.cc
+++ b/internal/ceres/schur_eliminator_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2019 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -284,10 +284,10 @@
BlockSparseMatrix matrix(bs);
double* values = matrix.mutable_values();
std::mt19937 prng;
- std::normal_distribution<double> standard_normal;
- for (int i = 0; i < matrix.num_nonzeros(); ++i) {
- values[i] = standard_normal(prng);
- }
+ std::normal_distribution<> standard_normal;
+ std::generate_n(values, matrix.num_nonzeros(), [&prng, &standard_normal] {
+ return standard_normal(prng);
+ });
Vector b(matrix.num_rows());
b.setRandom();
diff --git a/internal/ceres/sparse_cholesky_test.cc b/internal/ceres/sparse_cholesky_test.cc
index 73e1455..2397fdd 100644
--- a/internal/ceres/sparse_cholesky_test.cc
+++ b/internal/ceres/sparse_cholesky_test.cc
@@ -55,7 +55,8 @@
const int num_col_blocks,
const int min_col_block_size,
const int max_col_block_size,
- const double block_density) {
+ const double block_density,
+ std::mt19937& prng) {
// Create a random matrix
BlockSparseMatrix::RandomMatrixOptions options;
options.num_col_blocks = num_col_blocks;
@@ -66,7 +67,7 @@
options.min_row_block_size = 1;
options.max_row_block_size = max_col_block_size;
options.block_density = block_density;
- auto random_matrix = BlockSparseMatrix::CreateRandomMatrix(options);
+ auto random_matrix = BlockSparseMatrix::CreateRandomMatrix(options, prng);
// Add a diagonal block sparse matrix to make it full rank.
Vector diagonal = Vector::Ones(random_matrix->num_cols());
@@ -110,7 +111,8 @@
const int num_blocks,
const int min_block_size,
const int max_block_size,
- const double block_density) {
+ const double block_density,
+ std::mt19937& prng) {
LinearSolver::Options sparse_cholesky_options;
sparse_cholesky_options.sparse_linear_algebra_library_type =
sparse_linear_algebra_library_type;
@@ -120,7 +122,7 @@
sparse_cholesky->StorageType();
auto m = CreateRandomFullRankMatrix(
- num_blocks, min_block_size, max_block_size, block_density);
+ num_blocks, min_block_size, max_block_size, block_density, prng);
auto inner_product_computer = InnerProductComputer::Create(*m, storage_type);
inner_product_computer->Compute();
CompressedRowSparseMatrix* lhs = inner_product_computer->mutable_result();
@@ -184,7 +186,8 @@
num_blocks,
kMinBlockSize,
kMaxBlockSize,
- block_density);
+ block_density,
+ prng);
}
}
}
diff --git a/internal/ceres/subset_preconditioner_test.cc b/internal/ceres/subset_preconditioner_test.cc
index cd86695..5fb7f26 100644
--- a/internal/ceres/subset_preconditioner_test.cc
+++ b/internal/ceres/subset_preconditioner_test.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2017 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -31,6 +31,7 @@
#include "ceres/subset_preconditioner.h"
#include <memory>
+#include <random>
#include "Eigen/Dense"
#include "Eigen/SparseCore"
@@ -100,13 +101,13 @@
options.max_row_block_size = 4;
options.block_density = 0.9;
- m_ = BlockSparseMatrix::CreateRandomMatrix(options);
+ m_ = BlockSparseMatrix::CreateRandomMatrix(options, prng);
start_row_block_ = m_->block_structure()->rows.size();
// Ensure that the bottom part of the matrix has the same column
// block structure.
options.col_blocks = m_->block_structure()->cols;
- b_ = BlockSparseMatrix::CreateRandomMatrix(options);
+ b_ = BlockSparseMatrix::CreateRandomMatrix(options, prng);
m_->AppendRows(*b_);
// Create a Identity block diagonal matrix with the same column
@@ -132,6 +133,7 @@
std::unique_ptr<Preconditioner> preconditioner_;
Vector diagonal_;
int start_row_block_;
+ std::mt19937 prng;
};
TEST_P(SubsetPreconditionerTest, foo) {
diff --git a/internal/ceres/triplet_sparse_matrix.cc b/internal/ceres/triplet_sparse_matrix.cc
index 49d367a..37a779f 100644
--- a/internal/ceres/triplet_sparse_matrix.cc
+++ b/internal/ceres/triplet_sparse_matrix.cc
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -310,7 +310,8 @@
}
std::unique_ptr<TripletSparseMatrix> TripletSparseMatrix::CreateRandomMatrix(
- const TripletSparseMatrix::RandomMatrixOptions& options) {
+ const TripletSparseMatrix::RandomMatrixOptions& options,
+ std::mt19937& prng) {
CHECK_GT(options.num_rows, 0);
CHECK_GT(options.num_cols, 0);
CHECK_GT(options.density, 0.0);
@@ -319,7 +320,6 @@
std::vector<int> rows;
std::vector<int> cols;
std::vector<double> values;
- std::mt19937 prng;
std::uniform_real_distribution<double> uniform01(0.0, 1.0);
std::normal_distribution<double> standard_normal;
while (rows.empty()) {
diff --git a/internal/ceres/triplet_sparse_matrix.h b/internal/ceres/triplet_sparse_matrix.h
index c962455..549966a 100644
--- a/internal/ceres/triplet_sparse_matrix.h
+++ b/internal/ceres/triplet_sparse_matrix.h
@@ -1,5 +1,5 @@
// Ceres Solver - A fast non-linear least squares minimizer
-// Copyright 2015 Google Inc. All rights reserved.
+// Copyright 2022 Google Inc. All rights reserved.
// http://ceres-solver.org/
//
// Redistribution and use in source and binary forms, with or without
@@ -32,6 +32,7 @@
#define CERES_INTERNAL_TRIPLET_SPARSE_MATRIX_H_
#include <memory>
+#include <random>
#include <vector>
#include "ceres/crs_matrix.h"
@@ -135,7 +136,8 @@
// normally distributed and whose structure is determined by
// RandomMatrixOptions.
static std::unique_ptr<TripletSparseMatrix> CreateRandomMatrix(
- const TripletSparseMatrix::RandomMatrixOptions& options);
+ const TripletSparseMatrix::RandomMatrixOptions& options,
+ std::mt19937& prng);
// Load a triplet sparse matrix from a text file.
static std::unique_ptr<TripletSparseMatrix> CreateFromTextFile(FILE* file);
