Mantid
Loading...
Searching...
No Matches
DerivMinimizer.cpp
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
7//----------------------------------------------------------------------
8// Includes
9//----------------------------------------------------------------------
12
14
19double DerivMinimizer::fun(const gsl_vector *x, void *params) {
20 DerivMinimizer &minimizer = *static_cast<DerivMinimizer *>(params);
21 size_t n = minimizer.m_costFunction->nParams();
22 for (size_t i = 0; i < n; ++i) {
23 minimizer.m_costFunction->setParameter(i, gsl_vector_get(x, i));
24 }
25 std::shared_ptr<CostFunctions::CostFuncFitting> fitting =
26 std::dynamic_pointer_cast<CostFunctions::CostFuncFitting>(minimizer.m_costFunction);
27 if (fitting) {
28 fitting->getFittingFunction()->applyTies();
29 }
30 return minimizer.m_costFunction->val();
31}
32
38void DerivMinimizer::dfun(const gsl_vector *x, void *params, gsl_vector *g) {
39 DerivMinimizer &minimizer = *static_cast<DerivMinimizer *>(params);
40 size_t n = minimizer.m_costFunction->nParams();
41 for (size_t i = 0; i < n; ++i) {
42 minimizer.m_costFunction->setParameter(i, gsl_vector_get(x, i));
43 }
44 std::shared_ptr<CostFunctions::CostFuncFitting> fitting =
45 std::dynamic_pointer_cast<CostFunctions::CostFuncFitting>(minimizer.m_costFunction);
46 if (fitting) {
47 fitting->getFittingFunction()->applyTies();
48 }
49 std::vector<double> der(n);
50 minimizer.m_costFunction->deriv(der);
51 for (size_t i = 0; i < n; ++i) {
52 gsl_vector_set(g, i, der[i]);
53 }
54}
55
62void DerivMinimizer::fundfun(const gsl_vector *x, void *params, double *f, gsl_vector *g) {
63 DerivMinimizer &minimizer = *static_cast<DerivMinimizer *>(params);
64 size_t n = minimizer.m_costFunction->nParams();
65 for (size_t i = 0; i < n; ++i) {
66 minimizer.m_costFunction->setParameter(i, gsl_vector_get(x, i));
67 }
68 std::shared_ptr<CostFunctions::CostFuncFitting> fitting =
69 std::dynamic_pointer_cast<CostFunctions::CostFuncFitting>(minimizer.m_costFunction);
70 if (fitting) {
71 fitting->getFittingFunction()->applyTies();
72 }
73 std::vector<double> der(n);
74 *f = minimizer.m_costFunction->valAndDeriv(der);
75 for (size_t i = 0; i < n; ++i) {
76 gsl_vector_set(g, i, der[i]);
77 }
78}
79
82 : m_gslSolver(nullptr), m_x(nullptr), m_stopGradient(1e-3), m_stepSize(0.1), m_tolerance(0.0001) {
84}
85
91DerivMinimizer::DerivMinimizer(const double stepSize, const double tolerance)
92 : m_gslSolver(nullptr), m_x(nullptr), m_stopGradient(1e-3), m_stepSize(stepSize), m_tolerance(tolerance) {
94}
95
97 m_gslMultiminContainer.f = nullptr;
98 m_gslMultiminContainer.df = nullptr;
99 m_gslMultiminContainer.fdf = nullptr;
101 m_gslMultiminContainer.params = nullptr;
102}
103
108 if (m_gslSolver != nullptr) {
109 gsl_multimin_fdfminimizer_free(m_gslSolver);
110 gsl_vector_free(m_x);
111 }
112}
113
119void DerivMinimizer::initialize(API::ICostFunction_sptr function, size_t maxIterations) {
120 UNUSED_ARG(maxIterations);
121 m_costFunction = function;
126 m_gslMultiminContainer.params = this;
127
128 m_gslSolver = gsl_multimin_fdfminimizer_alloc(getGSLMinimizerType(), m_gslMultiminContainer.n);
129
130 size_t nParams = m_costFunction->nParams();
131 // Starting point
132 m_x = gsl_vector_alloc(nParams);
133 for (size_t i = 0; i < nParams; ++i) {
134 gsl_vector_set(m_x, i, m_costFunction->getParameter(i));
135 }
136
137 gsl_multimin_fdfminimizer_set(m_gslSolver, &m_gslMultiminContainer, m_x, m_stepSize, m_tolerance);
138}
139
144bool DerivMinimizer::iterate(size_t /*iteration*/) {
145 if (m_gslSolver == nullptr) {
146 throw std::runtime_error("Minimizer " + this->name() + " was not initialized.");
147 }
148 int status = gsl_multimin_fdfminimizer_iterate(m_gslSolver);
149 if (status) {
150 m_errorString = gsl_strerror(status);
151 return false;
152 }
153 status = gsl_multimin_test_gradient(m_gslSolver->gradient, m_stopGradient);
154 if (status != GSL_CONTINUE) {
155 m_errorString = gsl_strerror(status);
156 return false;
157 }
158 return true;
159}
160
166 if (value <= 0) {
167 throw std::invalid_argument("Gradient norm must be a positive number");
168 }
170}
171
173
174} // namespace Mantid::CurveFitting::FuncMinimisers
double value
The value of the point.
Definition: FitMW.cpp:51
double tolerance
#define UNUSED_ARG(x)
Function arguments are sometimes unused in certain implmentations but are required for documentation ...
Definition: System.h:64
std::string m_errorString
Error string.
virtual std::string name() const =0
Get name of minimizer.
A wrapper around the GSL functions implementing a minimizer using derivatives.
void setStopGradient(const double value)
Set maximum value of the gradient at which iterations can stop.
static void fundfun(const gsl_vector *x, void *params, double *f, gsl_vector *g)
Used by the GSL.
void initialize(API::ICostFunction_sptr function, size_t maxIterations=0) override
Initialize minimizer, i.e. pass a function to minimize.
gsl_multimin_function_fdf m_gslMultiminContainer
GSL container.
double m_stopGradient
the norm of the gradient at which iterations stop
bool iterate(size_t) override
Do one iteration.
static void dfun(const gsl_vector *x, void *params, gsl_vector *g)
Used by the GSL.
static double fun(const gsl_vector *x, void *params)
Used by the GSL.
gsl_multimin_fdfminimizer * m_gslSolver
pointer to the GSL solver doing the work
void initGSLMMin()
simply init the values for the gsl minimizer
virtual const gsl_multimin_fdfminimizer_type * getGSLMinimizerType()=0
Return a concrete type to initialize m_gslSolver with.
double costFunctionVal() override
Return current value of the cost function.
API::ICostFunction_sptr m_costFunction
Function to minimize.
gsl_vector * m_x
GSL vector with function parameters.
std::shared_ptr< ICostFunction > ICostFunction_sptr
define a shared pointer to a cost function
Definition: ICostFunction.h:60