GoogleGit

blob: 1c2e857989b4ce3abbefebfcadf0bda89d57371d [file] [log] [blame]
  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: sameeragarwal@google.com (Sameer Agarwal)
  30. //
  31. // An example program that minimizes Powell's singular function.
  32. //
  33. // F = 1/2 (f1^2 + f2^2 + f3^2 + f4^2)
  34. //
  35. // f1 = x1 + 10*x2;
  36. // f2 = sqrt(5) * (x3 - x4)
  37. // f3 = (x2 - 2*x3)^2
  38. // f4 = sqrt(10) * (x1 - x4)^2
  39. //
  40. // The starting values are x1 = 3, x2 = -1, x3 = 0, x4 = 1.
  41. // The minimum is 0 at (x1, x2, x3, x4) = 0.
  42. //
  43. // From: Testing Unconstrained Optimization Software by Jorge J. More, Burton S.
  44. // Garbow and Kenneth E. Hillstrom in ACM Transactions on Mathematical Software,
  45. // Vol 7(1), March 1981.
  46. #include <vector>
  47. #include "ceres/ceres.h"
  48. #include "gflags/gflags.h"
  49. #include "glog/logging.h"
  50. using ceres::AutoDiffCostFunction;
  51. using ceres::CostFunction;
  52. using ceres::Problem;
  53. using ceres::Solver;
  54. using ceres::Solve;
  55. struct F1 {
  56. template <typename T> bool operator()(const T* const x1,
  57. const T* const x2,
  58. T* residual) const {
  59. // f1 = x1 + 10 * x2;
  60. residual[0] = x1[0] + T(10.0) * x2[0];
  61. return true;
  62. }
  63. };
  64. struct F2 {
  65. template <typename T> bool operator()(const T* const x3,
  66. const T* const x4,
  67. T* residual) const {
  68. // f2 = sqrt(5) (x3 - x4)
  69. residual[0] = T(sqrt(5.0)) * (x3[0] - x4[0]);
  70. return true;
  71. }
  72. };
  73. struct F3 {
  74. template <typename T> bool operator()(const T* const x2,
  75. const T* const x4,
  76. T* residual) const {
  77. // f3 = (x2 - 2 x3)^2
  78. residual[0] = (x2[0] - T(2.0) * x4[0]) * (x2[0] - T(2.0) * x4[0]);
  79. return true;
  80. }
  81. };
  82. struct F4 {
  83. template <typename T> bool operator()(const T* const x1,
  84. const T* const x4,
  85. T* residual) const {
  86. // f4 = sqrt(10) (x1 - x4)^2
  87. residual[0] = T(sqrt(10.0)) * (x1[0] - x4[0]) * (x1[0] - x4[0]);
  88. return true;
  89. }
  90. };
  91. DEFINE_string(minimizer, "trust_region",
  92. "Minimizer type to use, choices are: line_search & trust_region");
  93. int main(int argc, char** argv) {
  94. CERES_GFLAGS_NAMESPACE::ParseCommandLineFlags(&argc, &argv, true);
  95. google::InitGoogleLogging(argv[0]);
  96. double x1 = 3.0;
  97. double x2 = -1.0;
  98. double x3 = 0.0;
  99. double x4 = 1.0;
  100. Problem problem;
  101. // Add residual terms to the problem using the using the autodiff
  102. // wrapper to get the derivatives automatically. The parameters, x1 through
  103. // x4, are modified in place.
  104. problem.AddResidualBlock(new AutoDiffCostFunction<F1, 1, 1, 1>(new F1),
  105. NULL,
  106. &x1, &x2);
  107. problem.AddResidualBlock(new AutoDiffCostFunction<F2, 1, 1, 1>(new F2),
  108. NULL,
  109. &x3, &x4);
  110. problem.AddResidualBlock(new AutoDiffCostFunction<F3, 1, 1, 1>(new F3),
  111. NULL,
  112. &x2, &x3);
  113. problem.AddResidualBlock(new AutoDiffCostFunction<F4, 1, 1, 1>(new F4),
  114. NULL,
  115. &x1, &x4);
  116. Solver::Options options;
  117. LOG_IF(FATAL, !ceres::StringToMinimizerType(FLAGS_minimizer,
  118. &options.minimizer_type))
  119. << "Invalid minimizer: " << FLAGS_minimizer
  120. << ", valid options are: trust_region and line_search.";
  121. options.max_num_iterations = 100;
  122. options.linear_solver_type = ceres::DENSE_QR;
  123. options.minimizer_progress_to_stdout = true;
  124. std::cout << "Initial x1 = " << x1
  125. << ", x2 = " << x2
  126. << ", x3 = " << x3
  127. << ", x4 = " << x4
  128. << "\n";
  129. // Run the solver!
  130. Solver::Summary summary;
  131. Solve(options, &problem, &summary);
  132. std::cout << summary.FullReport() << "\n";
  133. std::cout << "Final x1 = " << x1
  134. << ", x2 = " << x2
  135. << ", x3 = " << x3
  136. << ", x4 = " << x4
  137. << "\n";
  138. return 0;
  139. }