Re #4158. Debugged new Levenberg-Marquardt minimizer.

Code/Mantid/Framework/API/inc/MantidAPI/FunctionValues.h

@@ -77,6 +77,7 @@ class MANTID_API_DLL FunctionValues
   void setFitWeights(const double& value);
  /// get a fitting weight
   double getFitWeight(size_t i) const;
+  void setFitDataFromCalculated(const FunctionValues& values);
 protected:
   std::vector<double> m_calculated; ///< buffer for calculated values
   std::vector<double> m_data;    ///< buffer for fit data

Code/Mantid/Framework/API/src/FunctionValues.cpp

@@ -163,5 +163,13 @@ namespace API
     return m_weights[i];
   }
 
+  /**
+   * Set fitting data copied from other FunctionValues' calculated values.
+   */
+  void FunctionValues::setFitDataFromCalculated(const FunctionValues& values)
+  {
+    m_data.assign(values.m_calculated.begin(),values.m_calculated.end());
+  }
+
 } // namespace API
 } // namespace Mantid

Code/Mantid/Framework/CurveFitting/CMakeLists.txt

@@ -153,6 +153,7 @@ set ( TEST_FILES
 	test/GaussianTest.h
 #	test/IkedaCarpenterPVTest.h
 	test/LeastSquaresTest.h
+	test/LevenbergMarquardtTest.h
 	test/LinearBackgroundTest.h
 #	test/LinearTest.h
 #	test/LogNormalTest.h

Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/CostFuncFitting.h

@@ -7,8 +7,7 @@
 #include "MantidAPI/ICostFunction.h"
 #include "MantidAPI/IFunction.h"
 #include "MantidKernel/Matrix.h"
-
-//#include <gsl/gsl_matrix.h>
+#include "MantidCurveFitting/GSLVector.h"
 
 namespace Mantid
 {
@@ -58,6 +57,9 @@ class DLLExport CostFuncFitting : public API::ICostFunction
   virtual void setFittingFunction(API::IFunction_sptr function, 
     API::FunctionDomain_sptr domain, API::FunctionValues_sptr values);
 
+  /// Get fitting function.
+  virtual API::IFunction_sptr getFittingFunction(){return m_function;}
+
   /// Calculates covariance matrix
   /// @param covar :: Returned covariance matrix, here as 
   /// @param epsrel :: Is used to remove linear-dependent columns

Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/CostFuncLeastSquares.h

@@ -64,7 +64,7 @@ class DLLExport CostFuncLeastSquares : public CostFuncFitting
   /// @return :: The value of the function
   virtual double valAndDeriv(std::vector<double>& der) const;
 
-  virtual double valDerivHessian(GSLVector& der, GSLMatrix& hessian) const;
+  virtual double valDerivHessian(GSLVector& der, GSLMatrix& hessian, bool evalFunction = true) const;
 
 private:
   /// name of this minimizer

Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/GSLJacobian.h

@@ -1,5 +1,5 @@
-#ifndef MANTID_CURVEFITTING_GSLFUNCTIONS_H_
-#define MANTID_CURVEFITTING_GSLFUNCTIONS_H_
+#ifndef MANTID_CURVEFITTING_GSLJACOBIAN_H_
+#define MANTID_CURVEFITTING_GSLJACOBIAN_H_
 
 #include "MantidAPI/Jacobian.h"
 #include "MantidAPI/IFunction.h"
@@ -107,4 +107,4 @@ namespace Mantid
   } // namespace CurveFitting
 } // namespace Mantid
 
-#endif /*MANTID_CURVEFITTING_GSLFUNCTIONS_H_*/
+#endif /*MANTID_CURVEFITTING_GSLJACOBIAN_H_*/

Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/GSLMatrix.h

@@ -42,32 +42,47 @@ namespace Mantid
     gsl_matrix * m_matrix;
 
   public:
+    /// Constructor
+    GSLMatrix():m_matrix(NULL){}
     /// Constructor
     GSLMatrix(const size_t nx, const size_t ny)
     {
       m_matrix = gsl_matrix_alloc(nx,ny);
     }
 
+    /// Copy constructor
+    GSLMatrix(const GSLMatrix& M)
+    {
+      m_matrix = gsl_matrix_alloc(M.size1(),M.size2());
+      gsl_matrix_memcpy(m_matrix,M.gsl());
+    }
+
     /// Destructor.
     ~GSLMatrix()
     {
-      gsl_matrix_free(m_matrix);
+      if (m_matrix)
+      {
+        gsl_matrix_free(m_matrix);
+      }
     }
 
-    /// Get the pointer to the GSL's jacobian
+    /// Get the pointer to the GSL matrix
     gsl_matrix * gsl(){return m_matrix;}
 
+    /// Get the const pointer to the GSL matrix
+    const gsl_matrix * gsl() const {return m_matrix;}
+
     void resize(const size_t nx, const size_t ny)
     {
       gsl_matrix_free(m_matrix);
       m_matrix = gsl_matrix_alloc(nx,ny);
     }
 
     /// First size of the matrix
-    size_t size1() const {return m_matrix->size1;}
+    size_t size1() const {return m_matrix? m_matrix->size1 : 0;}
 
     /// Second size of the matrix
-    size_t size2() const {return m_matrix->size2;}
+    size_t size2() const {return m_matrix? m_matrix->size2 : 0;}
 
     /// set an element
     void set(size_t i, size_t j, double value)
@@ -79,11 +94,51 @@ namespace Mantid
       }
     }
     /// get an element
-    double get(size_t i, size_t j)
+    double get(size_t i, size_t j) const
     {
       if (i < m_matrix->size1 && j < m_matrix->size2) return gsl_matrix_get(m_matrix,i,j);
       throw std::out_of_range("GSLMatrix indices are out of range.");
     }
+
+    // Set this matrix to identity matrix
+    void identity()
+    {
+      gsl_matrix_set_identity( m_matrix );
+    }
+
+    // Set all elements to zero
+    void zero()
+    {
+      gsl_matrix_set_zero( m_matrix );
+    }
+
+    // add a matrix to this
+    GSLMatrix& operator+=(const GSLMatrix& M)
+    {
+      gsl_matrix_add( m_matrix, M.gsl() );
+      return *this;
+    }
+
+    // add a constant to this matrix
+    GSLMatrix& operator+=(const double& d)
+    {
+      gsl_matrix_add_constant( m_matrix, d );
+      return *this;
+    }
+
+    // subtract a matrix from this
+    GSLMatrix& operator-=(const GSLMatrix& M)
+    {
+      gsl_matrix_sub( m_matrix, M.gsl() );
+      return *this;
+    }
+
+    // multiply this matrix by a number
+    GSLMatrix& operator*=(const double& d)
+    {
+      gsl_matrix_scale( m_matrix, d );
+      return *this;
+    }
   };
 
   } // namespace CurveFitting

Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/GSLVector.h

@@ -40,6 +40,8 @@ namespace Mantid
     gsl_vector * m_vector;
 
   public:
+    /// Constructor
+    GSLVector():m_vector(NULL){}
     /// Constructor
     GSLVector(const size_t n)
     {
@@ -49,7 +51,10 @@ namespace Mantid
     /// Destructor.
     ~GSLVector()
     {
-      gsl_vector_free(m_vector);
+      if (m_vector)
+      {
+        gsl_vector_free(m_vector);
+      }
     }
 
     /// Get the pointer to the GSL vector
@@ -62,7 +67,7 @@ namespace Mantid
     }
 
     /// Size of the vector
-    size_t size() const {return m_vector->size;}
+    size_t size() const {return m_vector? m_vector->size : 0;}
 
     /// set an element
     void set(size_t i, double value)
@@ -74,7 +79,7 @@ namespace Mantid
       }
     }
     /// get an element
-    double get(size_t i)
+    double get(size_t i) const
     {
       if (i < m_vector->size) return gsl_vector_get(m_vector,i);
       throw std::out_of_range("GSLVector index is out of range.");

Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/IFuncMinimizer.h

@@ -59,7 +59,7 @@ class DLLExport IFuncMinimizer
 
   /// Perform iteration with minimizer and return info about how well this went
   /// using the GSL status integer system. See gsl_errno.h for details.
-  virtual bool minimize(size_t maxIterations = 500);
+  virtual bool minimize(size_t maxIterations = 1000);
 
   virtual std::string getError() const {return m_errorString;}
 

Code/Mantid/Framework/CurveFitting/inc/MantidCurveFitting/LevenbergMarquardtMinimizer.h

@@ -5,6 +5,8 @@
 // Includes
 //----------------------------------------------------------------------
 #include "MantidCurveFitting/IFuncMinimizer.h"
+#include "MantidCurveFitting/GSLVector.h"
+#include "MantidCurveFitting/GSLMatrix.h"
 
 namespace Mantid
 {
@@ -43,7 +45,7 @@ class DLLExport LevenbergMarquardtMinimizer : public IFuncMinimizer
 {
 public:
   /// Constructor
-  LevenbergMarquardtMinimizer():IFuncMinimizer(){}
+  LevenbergMarquardtMinimizer();
   /// Name of the minimizer.
   std::string name() const {return "Levenberg-Marquardt";}
 
@@ -55,8 +57,24 @@ class DLLExport LevenbergMarquardtMinimizer : public IFuncMinimizer
   virtual double costFunctionVal();
 
 private:
-
+  /// Pointer to the cost function. Must be the least squares.
   boost::shared_ptr<CostFuncLeastSquares> m_leastSquares;
+  /// Relative tolerance.
+  double m_relTol;
+  /// The tau parameter in the Levenberg-Marquardt method.
+  double m_tau;
+  /// The damping mu parameter in the Levenberg-Marquardt method.
+  double m_mu;
+  /// The nu parameter in the Levenberg-Marquardt method.
+  double m_nu;
+  /// The rho parameter in the Levenberg-Marquardt method.
+  double m_rho;
+  /// To keep function value
+  double m_F;
+  /// To keep first derivatives
+  GSLVector m_der;
+  /// To keep second derivatives
+  GSLMatrix m_hessian;
 	/// Static reference to the logger class
 	static Kernel::Logger& g_log;
 };

Code/Mantid/Framework/CurveFitting/src/CostFuncFitting.cpp

@@ -3,6 +3,7 @@
 //----------------------------------------------------------------------
 #include "MantidCurveFitting/CostFuncFitting.h"
 #include "MantidCurveFitting/GSLJacobian.h"
+#include "MantidAPI/IConstraint.h"
 
 #include <gsl/gsl_multifit_nlin.h>
 
@@ -59,6 +60,11 @@ void CostFuncFitting::setFittingFunction(API::IFunction_sptr function,
     {
       m_indexMap.push_back(i);
     }
+    API::IConstraint* c = m_function->getConstraint(i);
+    if (c)
+    {
+      c->setParamToSatisfyConstraint();
+    }
   }
 }