@@ -49,7 +49,7 @@ FGraphSolve::~FGraphSolve() = default;
4949
5050
5151uint_t FGraphSolve::solve (optimMethod method, uint_t maxIters,
52- matData_t lambda , matData_t solutionTolerance, bool verbose)
52+ matData_t lambdaParam , matData_t solutionTolerance, bool verbose)
5353{
5454 /* *
5555 * 2800 2D nodes on M3500
@@ -60,7 +60,7 @@ uint_t FGraphSolve::solve(optimMethod method, uint_t maxIters,
6060 * 1.3959 % update values,
6161 *
6262 */
63- lambda_ = lambda ;
63+ lambda_ = lambdaParam ;
6464 solutionTolerance_ = solutionTolerance;
6565 verbose_ = verbose;
6666 time_profiles_.reset ();
@@ -168,7 +168,7 @@ uint_t FGraphSolve::optimize_levenberg_marquardt(uint_t maxIters)
168168 // LM trust region as described in Bertsekas (p.105)
169169
170170 // 0) parameter initialization
171- lambda_ = 1e-5 ;
171+ // lambda_ = 1e-5;
172172 // sigma reference to the fidelity of the model at the proposed solution \in [0,1]
173173 matData_t sigma1 (0.25 ), sigma2 (0.8 );// 0 < sigma1 < sigma2 < 1
174174 matData_t beta1 (2.0 ), beta2 (0.25 ); // lambda updates multiplier values, beta1 > 1 > beta2 >0
@@ -206,7 +206,7 @@ uint_t FGraphSolve::optimize_levenberg_marquardt(uint_t maxIters)
206206 }
207207
208208 // 1.3) check for convergence
209- if (deltaChi2 < solutionTolerance_)
209+ if (deltaChi2/currentChi2 < solutionTolerance_) // in ratio
210210 {
211211 if (verbose_)
212212 std::cout << " \n FGraphSolve::optimize_levenberg_marquardt: Converged Successfully" << std::endl;
@@ -218,9 +218,14 @@ uint_t FGraphSolve::optimize_levenberg_marquardt(uint_t maxIters)
218218 // f = chi2(x_k) - chi2(x_k + dx)
219219 // chi2(x_k) - m_k(dx)
220220 // where m_k is the quadratized model = ||r||^2 - dx'*J' r + 0.5 dx'(J'J + lambda*D2)dx
221- modelFidelity = deltaChi2 / (dx_.dot (b_) - 0.5 *dx_.dot (L_* dx_));
221+ // modelFidelity = deltaChi2 / (dx_.dot(b_) - 0.5*dx_.dot(L_* dx_));
222+ // Update, according to H.B. Nielson, Damping Parameter In Marquardt’s Method, Technical Report IMM-REP-1999-05, Dept. of Mathematical Modeling, Technical University Denmark
223+ // (J'J ) dx = J'r , so this is substituted.
224+ modelFidelity = 2.0 * deltaChi2 / (dx_.dot (b_) - dx_.dot (lambda_ * diagL_* dx_));
222225 if (verbose_)
226+ {
223227 std::cout << " model fidelity = " << modelFidelity << " and m_k = " << dx_.dot (b_) << std::endl;
228+ }
224229
225230 // 3) update lambda
226231 if (modelFidelity < sigma1)
@@ -613,4 +618,4 @@ std::vector<bool> FGraphSolve::get_factors_robust_mask()
613618 }
614619
615620 return results;
616- }
621+ }
0 commit comments