55 const MnStrategy &strategy,
unsigned int maxfcn,
double edmval)
const
78 print.
Warn(
"No free parameters.");
83 print.
Error(
"Minimum seed invalid.");
88 print.
Error(
"Initial matrix not pos.def.");
94 std::vector<MinimumState>
result;
101 print.
Info(
"Start iterating until Edm is <", edmval,
"with call limit =", maxfcn);
106 int maxfcn_eff = maxfcn;
114 print.
Debug(ipass > 0 ?
"Continue" :
"Start",
"iterating...");
119 if (min.HasReachedCallLimit()) {
120 print.
Warn(
"FunctionMinimum is invalid, reached function call limit");
126 if (!min.IsValid()) {
127 print.
Warn(
"FunctionMinimum is invalid after second try");
133 edm =
result.back().Edm();
136 if ((strategy.
Strategy() >= 2) || (strategy.
Strategy() == 1 && min.Error().Dcovar() > 0.05)) {
138 print.
Debug(
"MnMigrad will verify convergence and Error matrix; dcov =", min.Error().Dcovar());
144 print.
Info(
"After Hessian");
149 print.
Warn(
"Invalid Hessian - exit the minimization");
156 print.
Debug(
"New Edm", edm,
"Requested", edmval);
161 if (edm >= machineLimit) {
164 print.
Info(
"Tolerance not sufficient, continue minimization; "
166 edm,
"Required", edmval);
168 print.
Warn(
"Reached machine accuracy limit; Edm", edm,
"is smaller than machine limit", machineLimit,
169 "while", edmval,
"was requested");
181 maxfcn_eff =
int(maxfcn * 1.3);
189 if (edm > 10 * edmval) {
191 print.
Warn(
"No convergence; Edm", edm,
"is above tolerance", 10 * edmval);
197 if (min.IsAboveMaxEdm())
198 print.
Info(
"Edm has been re-computed after Hesse; Edm", edm,
"is now within tolerance");
202 print.
Debug(
"Minimum found", min);
208 std::vector<MinimumState> &
result,
unsigned int maxfcn,
225 double edm = initialState.
Edm();
227 print.
Debug(
"Initial State:",
"\n Parameter:", initialState.
Vec(),
"\n Gradient:", initialState.
Gradient().
Vec(),
243 step = -1. *
s0.Error().InvHessian() *
s0.Gradient().Vec();
246 "\n Internal parameters",
s0.Vec(),
"\n Newton step", step);
250 print.
Debug(
"all derivatives are zero - return current status");
258 print.
Warn(
"Matrix not pos.def, gdel =", gdel,
"> 0");
262 step = -1. *
s0.Error().InvHessian() *
s0.Gradient().Vec();
269 print.
Warn(
"gdel =", gdel);
281 if (std::fabs(pp.
Y() -
s0.Fval()) <= std::fabs(
s0.Fval()) * prec.
Eps()) {
283 print.
Warn(
"No improvement in line search");
296 print.
Debug(
"Result after line search :",
"\n x =", pp.
X(),
"\n Old Fval =",
s0.Fval(),
297 "\n New Fval =", pp.
Y(),
"\n NFcalls =", fcn.
NumOfCalls());
305 if (std::isnan(edm)) {
306 print.
Warn(
"Edm is NaN; stop iterations");
312 print.
Warn(
"Matrix not pos.def., try to make pos.def.");
318 print.
Warn(
"Matrix still not pos.def.; stop iterations");
328 print.
Debug(
"Updated new point:",
"\n Parameter:",
p.Vec(),
"\n Gradient:",
g.Vec(),
329 "\n InvHessian:",
e.Matrix(),
"\n Edm:", edm);
340 edm *= (1. + 3. *
e.Dcovar());
342 print.
Debug(
"Dcovar =",
e.Dcovar(),
"\tCorrected edm =", edm);
344 }
while (edm > edmval && fcn.
NumOfCalls() < maxfcn);
350 if (!
result.back().IsValid())
354 print.
Warn(
"Call limit exceeded");
359 if (edm < 10 * edmval) {
360 print.
Info(
"Edm is close to limit - return current minimum");
362 }
else if (edm < std::fabs(prec.
Eps2() *
result.back().Fval())) {
363 print.
Warn(
"Edm is limited by Machine accuracy - return current minimum");
366 print.
Warn(
"Iterations finish without convergence; Edm", edm,
"Requested", edmval);
373 print.
Debug(
"Exiting successfully;",
"Ncalls", fcn.
NumOfCalls(),
"FCN",
result.back().Fval(),
"Edm", edm,
374 "Requested", edmval);
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void gc
const MnAlgebraicVector & Vec() const
class holding the full result of the minimization; both internal and external (MnUserParameterState) ...
interface class for gradient calculators
unsigned int size() const
void TraceIteration(int iter, const MinimumState &state) const
virtual MinimumError Update(const MinimumState &, const MinimumParameters &, const FunctionGradient &) const =0
MinimumError keeps the inv.
bool HasReachedCallLimit() const
const MnAlgebraicSymMatrix & InvHessian() const
const MnAlgebraicVector & Vec() const
const MinimumParameters & Parameters() const
const MnMachinePrecision & Precision() const
const MinimumState & State() const
MinimumState keeps the information (position, Gradient, 2nd deriv, etc) after one minimization step (...
const MinimumError & Error() const
const MnAlgebraicVector & Vec() const
const FunctionGradient & Gradient() const
Wrapper class to FCNBase interface used internally by Minuit.
unsigned int NumOfCalls() const
API class for calculating the numerical covariance matrix (== 2x Inverse Hessian == 2x Inverse 2nd de...
Implements a 1-dimensional minimization along a given direction (i.e.
Sets the relative floating point (double) arithmetic precision.
double Eps() const
eps returns the smallest possible number so that 1.+eps > 1.
double Eps2() const
eps2 returns 2*sqrt(eps)
double Y() const
Accessor to the y (second) coordinate.
double X() const
Accessor to the x (first) coordinate.
Force the covariance matrix to be positive defined by adding extra terms in the diagonal.
void Debug(const Ts &... args)
void Error(const Ts &... args)
void Info(const Ts &... args)
void Warn(const Ts &... args)
API class for defining four levels of strategies: low (0), medium (1), high (2), very high (>=3); act...
unsigned int Strategy() const
void SetHessianForcePosDef(unsigned int flag)
void AddResult(std::vector< MinimumState > &result, const MinimumState &state) const
FunctionMinimum Minimum(const MnFcn &, const GradientCalculator &, const MinimumSeed &, const MnStrategy &, unsigned int, double) const override
const VariableMetricEDMEstimator & Estimator() const
const MinimumErrorUpdator & ErrorUpdator() const
double Estimate(const FunctionGradient &, const MinimumError &) const
int iterate(rng_state_t *X)
double inner_product(const LAVector &, const LAVector &)
This file contains a specialised ROOT message handler to test for diagnostic in unit tests.