first minimiser implementation

configure.ac

@@ -12,7 +12,7 @@ AM_INIT_AUTOMAKE([-Wall -Werror])
 AM_SILENT_RULES([yes])
 AC_CONFIG_HEADERS([config.h])
 
-# Checks for programs.
+# Checks for programs
 AC_PROG_CXX
 AC_PROG_AWK
 AC_PROG_CC
@@ -26,7 +26,7 @@ m4_ifdef([AM_PROG_AR],[AM_PROG_AR])
 
 LT_INIT
 
-# Option to enable SSE random generator
+# Options
 AC_ARG_ENABLE([SSE],
 	[AS_HELP_STRING([--enable-SSE],
 		[compiles SSE version of ranlux random generator])],
@@ -35,6 +35,18 @@ AC_ARG_ENABLE([SSE],
 				[Define to 1 if your CPU support SSE instructions.])],
 	[]
 )
+AC_ARG_WITH([Minuit2],
+    [AS_HELP_STRING([--with-Minuit2=prefix],
+		[try this for a non-standard install prefix of the Minuit2 library])],
+    [AM_CFLAGS="$AM_CFLAGS -I$with_Minuit2/include"]
+    [AM_CXXFLAGS="$AM_CXXFLAGS -I$with_Minuit2/include"]
+	[AM_LDFLAGS="$AM_LDFLAGS -L$with_Minuit2/lib"],
+	[]
+)
+
+CFLAGS="$AM_CFLAGS $CFLAGS"
+CXXFLAGS="$AM_CXXFLAGS $CXXFLAGS"
+LDFLAGS="$AM_LDFLAGS $LDFLAGS"
 
 # Get compilers informations
 AX_COMPILER_VENDOR
@@ -58,6 +70,22 @@ AC_DEFINE_UNQUOTED([GXX_VERSION],["$GXX_VERSION"],
 
 # Checks for libraries.
 AC_CHECK_LIB([m],[cos],[],[AC_MSG_ERROR([libm library not found])])
+SAVED_LDFLAGS=$LDFLAGS
+LDFLAGS="$LDFLAGS -lMinuit2"
+AC_MSG_CHECKING([for ROOT::Minuit2::BasicMinimumError in -lMinuit2]);
+AC_LINK_IFELSE(
+	[AC_LANG_PROGRAM([#include <Minuit2/BasicMinimumError.h>],
+	[ROOT::Minuit2::BasicMinimumError dummy(0)])],
+	[LIBS="$LIBS -lMinuit2"]
+	[AC_DEFINE([HAVE_MINUIT2],
+				[1],
+				[Define to 1 if you have the `Minuit2' library (-lMinuit2).])]
+    [have_minuit=true]
+	[AC_MSG_RESULT([yes])],
+	[have_minuit=false]
+    [AC_MSG_RESULT([no])])
+AM_CONDITIONAL([HAVE_MINUIT], [test x$have_minuit = xtrue])
+LDFLAGS=$SAVED_LDFLAGS
 
 # Checks for header files.
 AC_HEADER_STDC

examples/Makefile.am

@@ -18,6 +18,7 @@ noinst_PROGRAMS =           \
     exCompiledDoubleFunction\
     exMat                   \
     exMathInterpreter       \
+    exMin                   \
     exPlot                  \
     exRand
 
@@ -29,6 +30,10 @@ exMat_SOURCES                    = exMat.cpp
 exMat_CFLAGS                     = -g -O2
 exMat_LDFLAGS                    = -L../latan/.libs -llatan
 
+exMin_SOURCES                    = exMin.cpp
+exMin_CFLAGS                     = -g -O2
+exMin_LDFLAGS                    = -L../latan/.libs -llatan
+
 exMathInterpreter_SOURCES        = exMathInterpreter.cpp
 exMathInterpreter_CFLAGS         = -g -O2
 exMathInterpreter_LDFLAGS        = -L../latan/.libs -llatan

examples/exMin.cpp

@@ -0,0 +1,29 @@
+#include <iostream>
+#include <latan/CompiledFunction.hpp>
+#include <latan/MinuitMinimizer.hpp>
+
+using namespace std;
+using namespace Latan;
+
+int main(int argc, char* argv[])
+{
+    string source;
+    
+    if (argc != 2)
+    {
+        cerr << "usage: " << argv[0] << " <function>" << endl;
+        
+        return EXIT_FAILURE;
+    }
+    source = argv[1];
+    
+    CompiledDoubleFunction f(1, source);
+    MinuitMinimizer minimizer;
+    double min;
+    
+    minimizer.setVerbosity(Minimizer::Verbosity::Debug);
+    min = minimizer(f)(0);
+    cout << "function minimum = " << min << endl;
+    
+    return EXIT_SUCCESS;
+}

latan/Makefile.am

@@ -38,6 +38,7 @@ liblatan_la_SOURCES =   \
     MathInterpreter.cpp \
     MathParser.ypp      \
     MathLexer.lpp       \
+    Minimizer.cpp       \
     Plot.cpp            \
     RandGen.cpp         \
     Sample.cpp          \
@@ -53,11 +54,18 @@ liblatan_la_HEADERS =   \
     Mat.hpp             \
     Math.hpp            \
     MathInterpreter.hpp \
+    Minimizer.hpp       \
     ParserState.hpp     \
     Plot.hpp            \
     RandGen.hpp         \
     Sample.hpp          \
     StatArray.hpp
+    
+if HAVE_MINUIT
+liblatan_la_SOURCES += MinuitMinimizer.cpp
+liblatan_la_HEADERS += MinuitMinimizer.hpp
+endif
+
 liblatan_la_CFLAGS = $(COM_CFLAGS)
 liblatan_la_CXXFLAGS = $(COM_CXXFLAGS)
 

latan/Minimizer.cpp

@@ -0,0 +1,58 @@
+/*
+ * Minimizer.cpp, part of LatAnalyze 3
+ *
+ * Copyright (C) 2013 - 2014 Antonin Portelli
+ *
+ * LatAnalyze 3 is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * LatAnalyze 3 is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with LatAnalyze 3.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#include <latan/Minimizer.hpp>
+#include <latan/includes.hpp>
+
+using namespace std;
+using namespace Latan;
+
+/******************************************************************************
+ *                       Minimizer implementation                             *
+ ******************************************************************************/
+// constructor /////////////////////////////////////////////////////////////////
+Minimizer::Minimizer(Verbosity verbosity)
+: verbosity_(verbosity)
+{}
+
+// access //////////////////////////////////////////////////////////////////////
+unsigned int Minimizer::getDim(void) const
+{
+    return static_cast<unsigned int>(x_.size());
+}
+
+DVec & Minimizer::getState(void)
+{
+    return x_;
+}
+
+Minimizer::Verbosity Minimizer::getVerbosity(void) const
+{
+    return verbosity_;
+}
+
+void Minimizer::setInit(const DVec &x0)
+{
+    x_ = x0;
+}
+
+void Minimizer::setVerbosity(const Verbosity verbosity)
+{
+    verbosity_ = verbosity;
+}

latan/Minimizer.hpp

@@ -0,0 +1,64 @@
+/*
+ * Minimizer.hpp, part of LatAnalyze 3
+ *
+ * Copyright (C) 2013 - 2014 Antonin Portelli
+ *
+ * LatAnalyze 3 is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * LatAnalyze 3 is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with LatAnalyze 3.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#ifndef Latan_Minimizer_hpp_
+#define Latan_Minimizer_hpp_
+
+#include <latan/Global.hpp>
+#include <latan/Function.hpp>
+#include <latan/Mat.hpp>
+
+BEGIN_NAMESPACE
+
+/******************************************************************************
+ *                        Abstract minimizer class                            *
+ ******************************************************************************/
+
+class Minimizer
+{
+public:
+    enum class Verbosity
+    {
+        Silent = 0,
+        Normal = 1,
+        Debug  = 2
+    };
+public:
+    // constructor
+    Minimizer(Verbosity verbosity = Verbosity::Silent);
+    // destructor
+    virtual ~Minimizer(void) = default;
+    // access
+    unsigned int getDim(void) const;
+    Verbosity    getVerbosity(void) const;
+    void         setInit(const DVec &x0);
+    void         setVerbosity(const Verbosity verbosity);
+    // minimization
+    virtual const DVec & operator()(const DoubleFunction &f) = 0;
+protected:
+    // access
+    DVec & getState(void);
+private:
+    DVec      x_;
+    Verbosity verbosity_;
+};
+
+END_NAMESPACE
+
+#endif // Latan_Minimizer_hpp_

latan/MinuitMinimizer.cpp

@@ -0,0 +1,129 @@
+/*
+ * MinuitMinimizer.cpp, part of LatAnalyze 3
+ *
+ * Copyright (C) 2013 - 2014 Antonin Portelli
+ *
+ * LatAnalyze 3 is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * LatAnalyze 3 is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with LatAnalyze 3.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#include <latan/MinuitMinimizer.hpp>
+#include <latan/includes.hpp>
+#include <Minuit2/FCNBase.h>
+#include <Minuit2/FunctionMinimum.h>
+#include <Minuit2/MnMigrad.h>
+#include <Minuit2/MnPrint.h>
+#include <Minuit2/MnPlot.h>
+#include <Minuit2/MnScan.h>
+#include <Minuit2/MnSimplex.h>
+#include <Minuit2/ScanMinimizer.h>
+#include <Minuit2/SimplexMinimizer.h>
+#include <Minuit2/VariableMetricMinimizer.h>
+
+using namespace std;
+using namespace ROOT;
+using namespace Minuit2;
+using namespace Latan;
+
+/******************************************************************************
+ *                    MinuitMinimizer implementation                          *
+ ******************************************************************************/
+// MinuitFunction constructor //////////////////////////////////////////////////
+MinuitMinimizer::MinuitFunction::MinuitFunction(const DoubleFunction &f)
+: f_(&f)
+{}
+
+// MinuitFunction minuit members ///////////////////////////////////////////////
+double MinuitMinimizer::MinuitFunction::operator()
+    (const vector<double> &x) const
+{
+    return f_->evaluate(x);
+}
+
+double MinuitMinimizer::MinuitFunction::Up(void) const
+{
+    return 1.;
+}
+
+// minimization ////////////////////////////////////////////////////////////////
+const DVec & MinuitMinimizer::operator()(const DoubleFunction &f)
+{
+    DVec &x = getState();
+    Verbosity verbosity = getVerbosity();
+    
+    if (f.getNArg() != x.size())
+    {
+        x.conservativeResize(f.getNArg());
+    }
+    
+    // set parameters
+    MinuitFunction minuitF(f);
+    MnUserParameters parameters;
+    
+    for (unsigned int i = 0; i < x.size(); ++i)
+    {
+        parameters.Add("x_" + strFrom(i), x(i), 0.1*fabs(x(i)));
+    }
+    
+    // pre-minimization
+    MnMigrad migrad0(minuitF, parameters, 0);
+    FunctionMinimum min = migrad0();
+    
+    if (verbosity >= Verbosity::Debug)
+    {
+        cout << "-- MINUIT pre-minimization -----------------------------";
+        cout << min;
+        cout << "--------------------------------------------------------";
+        cout << endl;
+    }
+    for (unsigned int i = 0; i < x.size(); ++i)
+    {
+        parameters.SetValue(i, min.UserParameters().Value(i));
+        parameters.SetError(i, min.UserParameters().Error(i));
+    }
+    
+    // minimization and output
+    MnMigrad  migrad2(minuitF, parameters, 2);
+    
+    min = migrad2();
+    for (unsigned int i = 0; i < x.size(); ++i)
+    {
+        x(i) = min.UserParameters().Value(i);
+    }
+    if (verbosity >= Verbosity::Normal)
+    {
+        cout << "-- MINUIT minimization ---------------------------------";
+        cout << min;
+        cout << "--------------------------------------------------------";
+        cout << endl;
+    }
+    if (verbosity >= Verbosity::Debug)
+    {
+        vector<pair<double, double>> scanRes;
+        MnPlot plot;
+        MnScan scanner(minuitF, min.UserParameters(), 2);
+        
+        cout << "-- MINUIT scan -----------------------------------------";
+        cout << endl;
+        for (unsigned int i = 0; i < x.size(); i++)
+        {
+            cout << "Parameter x_" << i << endl;
+            scanRes = scanner.Scan(i);
+            plot(scanRes);
+        }
+        cout << "--------------------------------------------------------";
+        cout << endl;
+    }
+    
+    return x;
+}

latan/MinuitMinimizer.hpp

@@ -0,0 +1,61 @@
+/*
+ * MinuitMinimizer.hpp, part of LatAnalyze 3
+ *
+ * Copyright (C) 2013 - 2014 Antonin Portelli
+ *
+ * LatAnalyze 3 is free software: you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation, either version 3 of the License, or
+ * (at your option) any later version.
+ *
+ * LatAnalyze 3 is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with LatAnalyze 3.  If not, see <http://www.gnu.org/licenses/>.
+ */
+
+#ifndef Latan_MinuitMinimizer_hpp_
+#define Latan_MinuitMinimizer_hpp_
+
+#include <latan/Global.hpp>
+#include <latan/Function.hpp>
+#include <latan/Minimizer.hpp>
+#include <Minuit2/FCNBase.h>
+
+BEGIN_NAMESPACE
+
+/******************************************************************************
+ *                             Minuit minimizer                               *
+ ******************************************************************************/
+
+class MinuitMinimizer: public Minimizer
+{
+private:
+    class MinuitFunction: public ROOT::Minuit2::FCNBase
+    {
+    public:
+        // constructor
+        explicit MinuitFunction(const DoubleFunction &f);
+        // destructor
+        virtual ~MinuitFunction(void) = default;
+        // minuit members
+        virtual double operator()(const std::vector<double> &x) const;
+        virtual double Up(void) const;
+    private:
+        const DoubleFunction *f_;
+    };
+public:
+    // constructors
+    using Minimizer::Minimizer;
+    // destructor
+    virtual ~MinuitMinimizer(void) = default;
+    // minimization
+    virtual const DVec & operator()(const DoubleFunction &f);
+};
+
+END_NAMESPACE
+
+#endif // Latan_MinuitMinimizer_hpp_