Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
GSLMultiMinimizer.h
Go to the documentation of this file.
1// @(#)root/mathmore:$Id$
2// Author: L. Moneta Tue Dec 19 14:09:15 2006
3
4/**********************************************************************
5 * *
6 * Copyright (c) 2006 LCG ROOT Math Team, CERN/PH-SFT *
7 * *
8 * This library is free software; you can redistribute it and/or *
9 * modify it under the terms of the GNU General Public License *
10 * as published by the Free Software Foundation; either version 2 *
11 * of the License, or (at your option) any later version. *
12 * *
13 * This library is distributed in the hope that it will be useful, *
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
16 * General Public License for more details. *
17 * *
18 * You should have received a copy of the GNU General Public License *
19 * along with this library (see file COPYING); if not, write *
20 * to the Free Software Foundation, Inc., 59 Temple Place, Suite *
21 * 330, Boston, MA 02111-1307 USA, or contact the author. *
22 * *
23 **********************************************************************/
24
25// Header file for class GSLMultiMinimizer
26
27#ifndef ROOT_Math_GSLMultiMinimizer
28#define ROOT_Math_GSLMultiMinimizer
29
30#include "gsl/gsl_vector.h"
31#include "gsl/gsl_multimin.h"
32#include "gsl/gsl_version.h"
34
35#include "Math/Error.h"
36
37#include "Math/IFunction.h"
38
39#include <cassert>
40#include <string>
41
42namespace ROOT {
43
44 namespace Math {
45
46
47/**
48 GSLMultiMinimizer class , for minimizing multi-dimensional function
49 using derivatives
50
51 @ingroup MultiMin
52
53*/
55
56public:
57
58 /**
59 Default constructor
60 */
62 fMinimizer(0),
63 fType(0),
64 fVec(0)
65 {
66 switch(type)
67 {
69 fType = gsl_multimin_fdfminimizer_conjugate_fr;
70 break;
72 fType = gsl_multimin_fdfminimizer_conjugate_pr;
73 break;
75 fType = gsl_multimin_fdfminimizer_vector_bfgs;
76 break;
78#if (GSL_MAJOR_VERSION > 1) || ((GSL_MAJOR_VERSION == 1) && (GSL_MINOR_VERSION >= 9))
79 // bfgs2 is available only for v>= 1.9
80 fType = gsl_multimin_fdfminimizer_vector_bfgs2;
81#else
82 MATH_INFO_MSG("GSLMultiMinimizer","minimizer BFSG2 does not exist with this GSL version , use BFGS");
83 fType = gsl_multimin_fdfminimizer_vector_bfgs;
84#endif
85 break;
87 fType = gsl_multimin_fdfminimizer_steepest_descent;
88 break;
89 default:
90 fType = gsl_multimin_fdfminimizer_conjugate_fr;
91 break;
92 }
93
94 }
95
96 /**
97 Destructor
98 */
100 if (fMinimizer != 0 ) gsl_multimin_fdfminimizer_free(fMinimizer);
101 // can free vector (is copied inside)
102 if (fVec != 0) gsl_vector_free(fVec);
103 }
104
105private:
106 // usually copying is non trivial, so we make this unaccessible
107
108 /**
109 Copy constructor
110 */
112
113 /**
114 Assignment operator
115 */
117 if (this == &rhs) return *this; // time saving self-test
118 return *this;
119 }
120
121public:
122
123 /**
124 set the function to be minimize the initial minimizer parameters,
125 step size and tolerance in the line search
126 */
127 int Set(const ROOT::Math::IMultiGradFunction & func, const double * x, double stepSize, double tol) {
128 // create function wrapper
129 fFunc.SetFunction(func);
130 // create minimizer object (free previous one if already existing)
131 unsigned int ndim = func.NDim();
132 CreateMinimizer( ndim );
133 // set initial values
134 if (fVec != 0) gsl_vector_free(fVec);
135 fVec = gsl_vector_alloc( ndim );
136 std::copy(x,x+ndim, fVec->data);
137 assert(fMinimizer != 0);
138 return gsl_multimin_fdfminimizer_set(fMinimizer, fFunc.GetFunc(), fVec, stepSize, tol);
139 }
140
141 /// create the minimizer from the type and size
142 void CreateMinimizer(unsigned int n) {
143 if (fMinimizer) gsl_multimin_fdfminimizer_free(fMinimizer);
144 fMinimizer = gsl_multimin_fdfminimizer_alloc(fType, n);
145 }
146
147 std::string Name() const {
148 if (fMinimizer == 0) return "undefined";
149 return std::string(gsl_multimin_fdfminimizer_name(fMinimizer) );
150 }
151
152 int Iterate() {
153 if (fMinimizer == 0) return -1;
154 return gsl_multimin_fdfminimizer_iterate(fMinimizer);
155 }
156
157 /// x values at the minimum
158 double * X() const {
159 if (fMinimizer == 0) return 0;
160 gsl_vector * x = gsl_multimin_fdfminimizer_x(fMinimizer);
161 return x->data;
162 }
163
164 /// function value at the minimum
165 double Minimum() const {
166 if (fMinimizer == 0) return 0;
167 return gsl_multimin_fdfminimizer_minimum(fMinimizer);
168 }
169
170 /// gradient value at the minimum
171 double * Gradient() const {
172 if (fMinimizer == 0) return 0;
173 gsl_vector * g = gsl_multimin_fdfminimizer_gradient(fMinimizer);
174 return g->data;
175 }
176
177 /// restart minimization from current point
178 int Restart() {
179 if (fMinimizer == 0) return -1;
180 return gsl_multimin_fdfminimizer_restart(fMinimizer);
181 }
182
183 /// test gradient (ask from minimizer gradient vector)
184 int TestGradient(double absTol) const {
185 if (fMinimizer == 0) return -1;
186 gsl_vector * g = gsl_multimin_fdfminimizer_gradient(fMinimizer);
187 return gsl_multimin_test_gradient( g, absTol);
188 }
189
190 /// test gradient (require a vector gradient)
191 int TestGradient(const double * g, double absTol) const {
192 if (fVec == 0 ) return -1;
193 unsigned int n = fVec->size;
194 if (n == 0 ) return -1;
195 std::copy(g,g+n, fVec->data);
196 return gsl_multimin_test_gradient( fVec, absTol);
197 }
198
199
200private:
201
202 gsl_multimin_fdfminimizer * fMinimizer;
204 const gsl_multimin_fdfminimizer_type * fType;
205 // cached vector to avoid re-allocating every time a new one
206 mutable gsl_vector * fVec;
207
208};
209
210 } // end namespace Math
211
212} // end namespace ROOT
213
214
215#endif /* ROOT_Math_GSLMultiMinimizer */
#define MATH_INFO_MSG(loc, str)
Pre-processor macro to report messages which can be configured to use ROOT error or simply an std::io...
Definition Error.h:77
#define g(i)
Definition RSha256.hxx:105
int type
Definition TGX11.cxx:121
Wrapper for a multi-dimensional function with derivatives used in GSL multidim minimization algorithm...
void SetFunction(const FuncType &f)
Fill gsl function structure from a C++ Function class.
GSLMultiMinimizer class , for minimizing multi-dimensional function using derivatives.
GSLMultiMinimizer & operator=(const GSLMultiMinimizer &rhs)
Assignment operator.
GSLMultiMinDerivFunctionWrapper fFunc
int Set(const ROOT::Math::IMultiGradFunction &func, const double *x, double stepSize, double tol)
set the function to be minimize the initial minimizer parameters, step size and tolerance in the line...
const gsl_multimin_fdfminimizer_type * fType
gsl_multimin_fdfminimizer * fMinimizer
double Minimum() const
function value at the minimum
double * Gradient() const
gradient value at the minimum
int Restart()
restart minimization from current point
int TestGradient(const double *g, double absTol) const
test gradient (require a vector gradient)
GSLMultiMinimizer(ROOT::Math::EGSLMinimizerType type)
Default constructor.
GSLMultiMinimizer(const GSLMultiMinimizer &)
Copy constructor.
double * X() const
x values at the minimum
void CreateMinimizer(unsigned int n)
create the minimizer from the type and size
int TestGradient(double absTol) const
test gradient (ask from minimizer gradient vector)
Interface (abstract class) for multi-dimensional functions providing a gradient calculation.
Definition IFunction.h:327
virtual unsigned int NDim() const=0
Retrieve the dimension of the function.
EGSLMinimizerType
enumeration specifying the types of GSL minimizers
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
Namespace for new Math classes and functions.
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...