Ceres作為一個優化算法庫，在許多領域中有著至關重要的作用，比如slam系統中的優化問題-集束調整BA，就可以通過Ceres去實現，官方文檔地址：http://ceres-solver.org/nnls_tutorial.html#bundle-adjustment

本文主要是解析ceres中的LM算法過程，參考代碼地址：

https://github.com/ceres-solver/ceres-solver/tree/master/internal/ceres

?

一、主要流程

先貼個圖，LM算法的大概流程如下

? ? ? ? ? ? ? ? ? ? ?

可以看到，LM算法的輸入為（1）雅可比矩陣J(x)；（2）殘差向量f(x)；（3）待優化變量初值x0；（4）控制參數等

LM算法要求解的問題為：

其中為殘差函數，它的導函數為,二階導函數的近似為

分為幾個步驟：

（1）初始化：首先計算系數矩陣A和殘差向量g,初始化參數

（2）while循環：如果達到收斂條件就停止迭代

? ? ? ? ? （3）求解 （A+miu*I)dx = -g

? ? ? ? ? ? ?(4)? xnew = x+dx

? ? ? ? ? ?（5）? ?,其中Fx,Fxnew是所有殘差的平方和

? ? ? ? ? ? ? ??

? ? ? ? ? ? ?(6)如果第五步結果大于零，表示這個迭代是有效的，可以接受

? ? ? ? ? ? ? ? ? ? ? ? ?然后更新x = xnew;Fx =Fnew ,A = J'*J, g = J'*f? ??

? ? ? ? ? ? ?(7)否則這個dx得到的結果是無效的，收縮搜索半徑，相當于增大

ceres對應代碼：https://github.com/ceres-solver/ceres-solver/blob/master/internal/ceres/trust_region_minimizer.cc

以及：https://github.com/ceres-solver/ceres-solver/blob/master/internal/ceres/levenberg_marquardt_strategy.cc

void TrustRegionMinimizer::Minimize(const Minimizer::Options& options,double* parameters,Solver::Summary* solver_summary) {start_time_in_secs_ = WallTimeInSeconds();iteration_start_time_in_secs_ = start_time_in_secs_; //初始化Init(options, parameters, solver_summary); //得到初始的雅可比矩陣以及殘差矩陣RETURN_IF_ERROR_AND_LOG(IterationZero());// Create the TrustRegionStepEvaluator. The construction needs to be// delayed to this point because we need the cost for the starting// point to initialize the step evaluator.step_evaluator_.reset(new TrustRegionStepEvaluator(x_cost_,options_.use_nonmonotonic_steps? options_.max_consecutive_nonmonotonic_steps: 0)); //while循環：是否達到迭代終止條件while (FinalizeIterationAndCheckIfMinimizerCanContinue()) {iteration_start_time_in_secs_ = WallTimeInSeconds();iteration_summary_ = IterationSummary();iteration_summary_.iteration =solver_summary->iterations.back().iteration + 1; //計算(J'*J+D'*D/radius)*dx = -J'*fRETURN_IF_ERROR_AND_LOG(ComputeTrustRegionStep());if (!iteration_summary_.step_is_valid) { //如果當前迭代不是有效的，則收縮搜索半徑RETURN_IF_ERROR_AND_LOG(HandleInvalidStep());continue;}if (options_.is_constrained &&options_.max_num_line_search_step_size_iterations > 0) {// Use a projected line search to enforce the bounds constraints// and improve the quality of the step.DoLineSearch(x_, gradient_, x_cost_, &delta_);} //計算Xnew = x+dxComputeCandidatePointAndEvaluateCost();DoInnerIterationsIfNeeded(); //是否||dx||太小了，太小了直接終止if (ParameterToleranceReached()) {return;} //是否||Fx-Fxnew||太小了，太小了意味著殘差平方和更新不動直接終止if (FunctionToleranceReached()) {return;} //如果是有效的if (IsStepSuccessful()) { //更新x,雅可比矩陣，殘差向量等RETURN_IF_ERROR_AND_LOG(HandleSuccessfulStep());continue;} //否則收縮搜索半徑HandleUnsuccessfulStep();} }

上面為什么是信賴域算法過程，Ceres里面把LM算法和dogleg算法（也叫狗腿算法）集成到統一的框架下--信賴域算法框架，不同的是LM算法求解dx的過程和狗腿算法不同，下面是LM算法求解dx的過程以及搜索半徑的更新

TrustRegionStrategy::Summary LevenbergMarquardtStrategy::ComputeStep(const TrustRegionStrategy::PerSolveOptions& per_solve_options,SparseMatrix* jacobian,const double* residuals,double* step) {CHECK(jacobian != nullptr);CHECK(residuals != nullptr);CHECK(step != nullptr);//計算對角矩陣Dconst int num_parameters = jacobian->num_cols();if (!reuse_diagonal_) {if (diagonal_.rows() != num_parameters) {diagonal_.resize(num_parameters, 1);}jacobian->SquaredColumnNorm(diagonal_.data());for (int i = 0; i < num_parameters; ++i) {diagonal_[i] = std::min(std::max(diagonal_[i], min_diagonal_),max_diagonal_);}} // D = diag{J'*J}/radiuslm_diagonal_ = (diagonal_ / radius_).array().sqrt();LinearSolver::PerSolveOptions solve_options;solve_options.D = lm_diagonal_.data();solve_options.q_tolerance = per_solve_options.eta;solve_options.r_tolerance = -1.0;InvalidateArray(num_parameters, step); //求解 （A+D'*D)dx = -gLinearSolver::Summary linear_solver_summary =linear_solver_->Solve(jacobian, residuals, solve_options, step);if (linear_solver_summary.termination_type == LINEAR_SOLVER_FATAL_ERROR) {LOG(WARNING) << "Linear solver fatal error: "<< linear_solver_summary.message;} else if (linear_solver_summary.termination_type == LINEAR_SOLVER_FAILURE) {LOG(WARNING) << "Linear solver failure. Failed to compute a step: "<< linear_solver_summary.message;} else if (!IsArrayValid(num_parameters, step)) {LOG(WARNING) << "Linear solver failure. Failed to compute a finite step.";linear_solver_summary.termination_type = LINEAR_SOLVER_FAILURE;} else {VectorRef(step, num_parameters) *= -1.0;}reuse_diagonal_ = true;if (per_solve_options.dump_format_type == CONSOLE ||(per_solve_options.dump_format_type != CONSOLE &&!per_solve_options.dump_filename_base.empty())) {if (!DumpLinearLeastSquaresProblem(per_solve_options.dump_filename_base,per_solve_options.dump_format_type,jacobian,solve_options.D,residuals,step,0)) {LOG(ERROR) << "Unable to dump trust region problem."<< " Filename base: " << per_solve_options.dump_filename_base;}}TrustRegionStrategy::Summary summary;summary.residual_norm = linear_solver_summary.residual_norm;summary.num_iterations = linear_solver_summary.num_iterations;summary.termination_type = linear_solver_summary.termination_type;return summary; }//如果當前迭代是有效的，更新搜索半徑等參數 void LevenbergMarquardtStrategy::StepAccepted(double step_quality) {CHECK_GT(step_quality, 0.0);radius_ = radius_ / std::max(1.0 / 3.0,1.0 - pow(2.0 * step_quality - 1.0, 3));radius_ = std::min(max_radius_, radius_);decrease_factor_ = 2.0;reuse_diagonal_ = false; }//如果當前迭代是無效的，收縮搜索半徑等參數 void LevenbergMarquardtStrategy::StepRejected(double step_quality) {radius_ = radius_ / decrease_factor_;decrease_factor_ *= 2.0;reuse_diagonal_ = true; }

?

總結

以上是生活随笔為你收集整理的ceres之LM算法的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。