Merge

external/quadprog++/Array.cc

+// This file is part of QuadProg++: a C++ library implementing
+// the algorithm of Goldfarb and Idnani for the solution of a (convex) 
+// Quadratic Programming problem by means of an active-set dual method.
+// Copyright (C) 2007-2009 Luca Di Gaspero.
+// Copyright (C) 2009 Eric Moyer.  
+//
+// QuadProg++ is free software: you can redistribute it and/or modify
+// it under the terms of the GNU Lesser General Public License as published by
+// the Free Software Foundation, either version 3 of the License, or
+// (at your option) any later version.
+//
+// QuadProg++ 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 Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public License
+// along with QuadProg++. If not, see <http://www.gnu.org/licenses/>.
+
+#include "Array.hh"
+
+/**
+  Index utilities
+ */
+namespace QuadProgPP{
+
+std::set<unsigned int> seq(unsigned int s, unsigned int e)
+{
+	std::set<unsigned int> tmp;
+	for (unsigned int i = s; i <= e; i++)
+		tmp.insert(i);
+	
+	return tmp;
+}
+
+std::set<unsigned int> singleton(unsigned int i)
+{
+	std::set<unsigned int> tmp;
+	tmp.insert(i);
+	
+	return tmp;
+}
+
+}

external/quadprog++/QuadProg++.hh

+/*
+ 
+ The quadprog_solve() function implements the algorithm of Goldfarb and Idnani 
+ for the solution of a (convex) Quadratic Programming problem
+ by means of an active-set dual method.
+	 
+The problem is in the form:
+
+min 0.5 * x G x + g0 x
+s.t.
+    CE^T x + ce0 = 0
+    CI^T x + ci0 >= 0
+	 
+ The matrix and vectors dimensions are as follows:
+     G: n * n
+		g0: n
+				
+		CE: n * p
+	 ce0: p
+				
+	  CI: n * m
+   ci0: m
+
+     x: n
+ 
+ The function will return the cost of the solution written in the x vector or
+ std::numeric_limits::infinity() if the problem is infeasible. In the latter case
+ the value of the x vector is not correct.
+ 
+ References: D. Goldfarb, A. Idnani. A numerically stable dual method for solving
+             strictly convex quadratic programs. Mathematical Programming 27 (1983) pp. 1-33.
+
+ Notes:
+  1. pay attention in setting up the vectors ce0 and ci0. 
+	   If the constraints of your problem are specified in the form 
+	   A^T x = b and C^T x >= d, then you should set ce0 = -b and ci0 = -d.  
+  2. The matrix G is modified within the function since it is used to compute
+     the G = L^T L cholesky factorization for further computations inside the function. 
+     If you need the original matrix G you should make a copy of it and pass the copy
+     to the function.
+    
+ Author: Luca Di Gaspero
+  			 DIEGM - University of Udine, Italy
+				 l.digaspero@uniud.it
+				 http://www.diegm.uniud.it/digaspero/
+ 
+ The author will be grateful if the researchers using this software will
+ acknowledge the contribution of this function in their research papers.
+
+LICENSE
+
+This file is part of QuadProg++: a C++ library implementing
+the algorithm of Goldfarb and Idnani for the solution of a (convex) 
+Quadratic Programming problem by means of an active-set dual method.
+Copyright (C) 2007-2009 Luca Di Gaspero.
+Copyright (C) 2009 Eric Moyer.  
+
+QuadProg++ is free software: you can redistribute it and/or modify
+it under the terms of the GNU Lesser General Public License as published by
+the Free Software Foundation, either version 3 of the License, or
+(at your option) any later version.
+
+QuadProg++ 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 Lesser General Public License for more details.
+
+You should have received a copy of the GNU Lesser General Public License
+along with QuadProg++. If not, see <http://www.gnu.org/licenses/>.
+
+*/
+
+
+#ifndef _QUADPROGPP
+#define _QUADPROGPP
+
+#include "Array.hh"
+
+namespace QuadProgPP{
+  double solve_quadprog(Matrix<double>& G, Vector<double>& g0, 
+			const Matrix<double>& CE, const Vector<double>& ce0,  
+			const Matrix<double>& CI, const Vector<double>& ci0, 
+			Vector<double>& x);
+}
+
+#endif // #define _QUADPROGPP

external/quadprog++/quadprog++.pro

+TARGET          = quadprog++
+TEMPLATE        = lib
+CONFIG         *= static  \
+				  qt      \
+                  plugin  \
+                  eigen
+
+DESTDIR         = ../../build
+ 
+INCLUDEPATH    += ../.. 
+
+HEADERS         = QuadProg++.hh \
+                  Array.hh
+SOURCES         = QuadProg++.cc \
+                  Array.cc
+

sources/core/args.h

@@ -23,7 +23,7 @@ class arguments
 		// Constructor and destructor
 		arguments()
 		{
-		} ;
+        }
 		arguments(int argc, char** const argv)
 		{
 			std::string key ;
@@ -53,8 +53,10 @@ class arguments
 				}
 				_map.insert(std::pair<std::string, std::string>(key, data)) ;
 			}
-		} ;
-		~arguments() { } ;
+        }
+        ~arguments()
+        {
+        }
 
 		//! \brief is the elements in the command line ?
 		bool is_defined(const std::string& key) const
@@ -67,7 +69,7 @@ class arguments
 			{
 				return false ;
 			}
-		} ;
+        }
 		//! \brief access the element stored value
 		std::string operator[](const std::string& key) const
 		{
@@ -80,7 +82,7 @@ class arguments
                 //std::cerr << "Underfined request to key : \"" << key << "\"" << std::endl ;
 				return std::string() ;
 			}
-		} ;
+        }
 		//! \brief acces to the float value associated with the parameter
 		//! \a key.
 		//!
@@ -113,10 +115,13 @@ class arguments
 		vec get_vec(const std::string& key, int size, float default_value = 0.0f) const
 		{
 			vec res(size);
+            for(int i=0; i<size; ++i)
+                res[i] = default_value;
+
 			if(_map.count(key) > 0)
 			{
 				std::string s = _map.find(key)->second;
-				if(s[0] == '\[') // Is an array of type [a, b, c]
+				if(s[0] == '[') // Is an array of type [a, b, c]
 				{
 					int i = 0;
 					size_t pos = 1;

sources/core/clustering.cpp

+#include "clustering.h"
+
+#include <iostream>
+#include <limits>
+#include <string>
+#include <fstream>
+
+#ifdef OLD
+clustering::clustering(const data* d, const arguments& args)
+{
+    // List of the indices of the clustering
+    std::vector<int> indices;
+
+    vec findices = args.get_vec("cluster-dim", d->dimX(), -1.0f);
+    for(int i=0; i<findices.size() && findices[i] != -1.0f; ++i)
+    {
+        int id = (int)findices[i];
+        indices.push_back(id);
+
+        assert(id >= 0 && id < d->dimX());
+    }
+    std::cout << "<<DEBUG>> the resulting data will have " << indices.size() << " dimensions" << std::endl;
+
+    // We cannot have null cluster size of full rank.
+    assert(indices.size() > 0 && indices.size() < d->dimX());
+
+    // Fit the slice into the domain of definition of the data
+    vec compressed_vec  = args.get_vec("cluster-slice", d->dimX(), std::numeric_limits<float>::max());
+    for(int i=0; i<d->dimX(); ++i)
+    {
+        compressed_vec[i] = std::min<double>(compressed_vec[i], d->max()[i]);
+        compressed_vec[i] = std::max<double>(compressed_vec[i], d->min()[i]);
+    }
+
+    // Set the limit of the new domain and dimensions
+    _nX = indices.size();
+    _nY = d->dimY();
+    _min = vec(_nX);
+    _max = vec(_nX);
+    for(int i=0; i<_nX; ++i)
+    {
+        _min[i] = d->min()[indices[i]];
+        _max[i] = d->max()[indices[i]];
+    }
+
+    for(int i=0; i<d->size(); ++i)
+    {
+        vec p = d->get(i);
+        bool reject = false;
+        for(int i=0; i<d->dimX(); ++i)
+        {
+            if(is_in<int>(indices, i) == -1 && p[i] != compressed_vec[i])
+            {
+                reject = true;
+            }
+        }
+        if(reject)
+        {
+            continue;
+        }
+
+        vec e(dimX()+dimY());
+        for(int i=0; i<_nX; ++i)
+        {
+            e[i] = p[indices[i]];
+        }
+        for(int i=0; i<_nY; ++i)
+        {
+            e[_nX + i] = p[d->dimX() + i];
+        }
+        _data.push_back(e);
+    }
+    std::cout << "<<INFO>> clustering left " << _data.size() << " elements" << std::endl;
+
+    save("cluster.txt");
+}
+
+void clustering::save(const std::string& filename)
+{
+    std::ofstream file(filename.c_str(), std::ios_base::trunc);
+    file << "#DIM " << _nX << " " << _nY << std::endl;
+    for(int i=0; i<size(); ++i)
+    {
+        for(int j=0; j<_nX+_nY; ++j)
+        {
+            file << _data[i][j] << "\t";
+        }
+        file << std::endl;
+    }
+    file.close();
+}
+
+vec clustering::get(int i) const
+{
+    return _data[i];
+}
+
+vec clustering::operator[](int i) const
+{
+    return _data[i];
+}
+
+int clustering::size() const
+{
+    return _data.size();
+}
+
+vec clustering::min() const
+{
+    return _min;
+}
+
+vec clustering::max() const
+{
+    return _max;
+}
+#else
+template<class T> void clustering(const T* in_data, int _nY, params::input in_param, params::input out_param, std::vector<vec>& out_data)
+{
+
+    // Set the dimensions
+    const int in_nX  = params::dimension(in_param);
+    const int out_nX = params::dimension(out_param);
+
+#ifdef DEBUG
+    std::cout << " in dim = " <<  in_nX << std::endl;
+    std::cout << "out dim = " << out_nX << std::endl;
+#endif
+
+    for(int i=0; i<in_data->size(); ++i)
+    {
+        vec p = in_data->get(i);
+        vec e (out_nX + _nY);
+
+        // Fill the input part of the vector
+        params::convert(&p[0], in_param, out_param, &e[0]);
+
+        // Search for duplicates
+		  
+        bool already_exist = false;
+        for(unsigned int j=0; j<out_data.size(); ++j)
+        {
+            vec test_e = out_data[j];
+
+            double dist = 0.0;
+            for(int k=0; k<out_nX; ++k)
+            {
+                const double temp = test_e[k]-e[k];
+                dist += temp*temp;
+            }
+
+            if(dist < 1.0E-5)
+                already_exist = true;
+        }
+
+        if(!already_exist)
+        {
+            // Fill the output part of the vector
+            for(int j=0; j<_nY; ++j)
+                e[out_nX + j] = p[in_nX + j];
+
+#ifdef DEBUG
+            std::cout << " in = " << p << std::endl;
+            std::cout << "out = " << e << std::endl;
+#endif
+            out_data.push_back(e);
+        }
+    }
+
+}
+#endif

sources/core/clustering.h

+#pragma once
+
+#include "common.h"
+#include "params.h"
+
+#include <vector>
+
+#ifdef OLD
+class clustering : public data
+{
+    public: // methods
+
+        //! \brief constructor loading a full dimension data and clustering
+        //! it into a low dimension one.
+        clustering::clustering(const data* d, const arguments& args);
+
+        //! \brief the clustering class can save a clustering to a file.
+        virtual void save(const std::string& filename);
+
+        //! \brief the clustering class cannot load from a file. It requires
+        //! a complete object to be created.
+        virtual void load(const std::string& filename)
+        {
+            throw;
+        }
+        //! \brief the clustering class cannot load from a file. It requires
+        //! a complete object to be created.
+        virtual void load(const std::string& filename, const arguments& args)
+        {
+            throw;
+        }
+
+        //! \brief aces to data in linear order
+        virtual vec get(int i) const ;
+        //! \brief aces to data in linear order
+        virtual vec operator[](int i) const ;
+
+        //! \brief return the size of the data after clusterization
+        int size() const;
+
+        //! \brief get min input space values
+        virtual vec min() const ;
+        //! \brief get max input space values
+        virtual vec max() const ;
+
+    protected:
+        std::vector<vec> _data;
+        vec _min, _max;
+};
+#else
+    template<class T> void clustering(const T* in_data , int nY, params::input in_param, params::input out_param, std::vector<vec>& out_data);
+#endif
+
+#include "clustering.cpp"

sources/core/common.h

@@ -4,6 +4,7 @@
 #include <iostream>
 #include <cassert>
 #include <cmath>
+#include <cstring>
 
 /*! \brief A core implementation of a vector of double.
  *  \ingroup core
@@ -28,7 +29,15 @@ class vec : public std::vector<double>
 		}
 		virtual ~vec() 
 		{
-		}
+        }
+
+        //! \brief get a subpart of the vector
+        vec subvector(int start, int n) const
+        {
+            vec res(n);
+            memcpy(&res[0], &this->at(start), n*sizeof(double));
+            return res;
+        }
 
 		//! \brief copy operator. It resize the left operand to the size of the 
 		//! right operand.
@@ -41,7 +50,7 @@ class vec : public std::vector<double>
 			}
 			return *this ;
 		}
-		
+
 		// Mathematical operators
 		//
 		friend vec operator-(const vec& a)
@@ -137,6 +146,16 @@ class vec : public std::vector<double>
 			}
 			return c ;
 		}
+        friend double norm(const vec& a)
+        {
+            double norm = 0.0 ;
+            for(unsigned int i=0; i<a.size(); ++i)
+            {
+                norm += a[i]*a[i];
+            }
+            return sqrt(norm);
+        }
+
 		friend vec normalize(const vec& a)
 		{
 			vec b(a.size());
@@ -169,14 +188,14 @@ class vec : public std::vector<double>
         friend bool operator<(const vec& a, const vec& b)
         {
             bool lessthan = true ;
-            for(int i=0; i<a.size(); ++i)
+            for(unsigned int i=0; i<a.size(); ++i)
                 lessthan &= a[i] < b[i];
             return lessthan;
         }
         friend bool operator>(const vec& a, const vec& b)
         {
             bool greatthan = true ;
-            for(int i=0; i<a.size(); ++i)
+            for(unsigned int i=0; i<a.size(); ++i)
                 greatthan &= a[i] > b[i];
             return greatthan;
         }
@@ -200,7 +219,36 @@ class vec : public std::vector<double>
 			out << "]" ;
 
 			return out ;
-		} ;
+        }
 
 
 } ;
+
+//! \brief locate the first index of value v in vector vec. Complexity in
+//! O(n) is the worst case.
+template<typename T> int is_in(std::vector<T> ve, T v)
+{
+    int res = -1;
+    for(unsigned int i=0; i<ve.size(); ++i)
+    {
+        if(ve[i] == v)
+            return i;
+    }
+
+    return res;
+}
+
+#define NOT_IMPLEMENTED() \
+std::cerr << "<<ERROR>> not implemented " << __FILE__ \
+          << ":" << __LINE__ << std::endl; \
+throw
+
+#ifdef WIN32
+#define M_PI 3.14159265
+#endif
+
+#ifdef WIN32
+#define ALTA_DLL_EXPORT extern "C" __declspec(dllexport)
+#else
+#define ALTA_DLL_EXPORT extern "C"
+#endif
\ No newline at end of file

sources/core/core.pro

@@ -4,6 +4,8 @@ CONFIG  *= static \
 
 DESTDIR  = ../build
 
+QMAKE_CXXFLAGS += -std=c++0x
+
 HEADERS  = args.h               \
            common.h             \
 		 	  data.h               \
@@ -12,8 +14,12 @@ HEADERS  = args.h               \
  			  plugins_manager.h    \
 			  vertical_segment.h   \
 			  rational_function.h  \
-			  params.h
+			  params.h             \
+           clustering.h
 
-SOURCES  = plugins_manager.cpp  \
-           vertical_segment.cpp \
-			  rational_function.cpp
+SOURCES  = plugins_manager.cpp   \
+           vertical_segment.cpp  \
+		   rational_function.cpp \
+		   params.cpp			 \
+           function.cpp          \
+#          clustering.cpp

sources/core/data.h

@@ -2,53 +2,104 @@
 
 #include <string>
 #include <utility>
+#include <iostream>
+#include <limits>
+#include <fstream>
 
 #include <QtPlugin>
 
 #include "common.h"
 #include "args.h"
 #include "params.h"
+#include "clustering.h"
 
 /*! \brief A data object. Allows to load data from files.
  *  \ingroup core
  */
 class data
 {
-		public: // methods
+	public: // methods
 
 		// Load data from a file
 		virtual void load(const std::string& filename) = 0 ;
 		virtual void load(const std::string& filename, const arguments& args) = 0 ;
 
+        // Save the data to a file
+        virtual void save(const std::string& filename) const
+        {
+            std::ofstream file(filename.c_str(), std::ios_base::trunc);
+            file << "#DIM " << _nX << " " << _nY << std::endl;
+            for(int i=0; i<size(); ++i)
+            {
+                vec x = this->get(i);
+                for(int j=0; j<_nX+_nY; ++j)
+                {
+                    file << x[j] << "\t";
+                }
+                file << std::endl;