NelderMeadTestFunctions.hpp

// Unfit: Data fitting and optimization software
//
// Copyright (C) 2012- Dr Martin Buist & Dr Alberto Corrias
// Contacts: martin.buist _at_ nus.edu.sg; alberto _at_ nus.edu.sg
//
// See the 'Contributors' file for a list of those who have contributed
// to this work.
//
// This program 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.
//
// This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.
//
#ifndef UNFIT_UNITTESTS_NELDERMEADTESTFUNCTIONS_HPP_
#define UNFIT_UNITTESTS_NELDERMEADTESTFUNCTIONS_HPP_

#include <cmath>
#include <limits>
#include <vector>
#include "GenericCostFunction.hpp"

namespace Unfit
{
namespace UnitTests
{
class LineToInfinity : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    if (x[0] > 3.0) {
        ret[0] = sqrt(fabs(x[0]));
    }
    else {
        ret[0] = std::numeric_limits<double>::infinity();
    }
    return ret;
  }
};

class InfinityInBetween : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    if (x[0] > 1.0) {
        ret[0] = x[0];
    }
    else if (x[0] < -1.0) {
        ret[0] = x[0];
    }
    else {
        ret[0] = std::numeric_limits<double>::infinity();
    }
    return ret;
  }
};

class Infinity2D : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> f_i(3, 0.0);
    if ((x[0] < 1.26 && x[0] > 1.24) && (x[1] < 0.51 && x[1] > 0.49)) {
      f_i[0] = std::numeric_limits<double>::infinity();
      f_i[1] = std::numeric_limits<double>::infinity();
      f_i[2] = std::numeric_limits<double>::infinity();
    }
    else if (x[0] < 0.9 || x[0] > 2.1 || x[1] < 0.9 || x[1] > 2.1) {
      f_i[0] = 3.5;
      f_i[1] = 3.5;
      f_i[2] = 3.5;
    }
    else {
      f_i[0] = x[0] + x[1];
      f_i[1] = x[0] + x[1];
      f_i[2] = x[0] + x[1];
    }
    return f_i;
  }
};

class SampleCostFunction1 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = sqrt(fabs(x[0]*x[0] - 4*x[0] + x[1]*x[1] - x[1] - x[0]*x[1]));
    return ret;
  }
};

// /**
// * The Sample function 2 is defined as
// *
// *   (y^2 + 4xy +x^3 - 2x).
// *
// * Number of dimensions = 2
// * The global minimum is: (2.896805253, -5.793610507)
// * Initial guess: (0, 0.8)
// *
// * Reference: http://www.wolframalpha.com/input/
// * ?i=%28y^2+%2B+4xy+%2Bx^3+-+2x%29
// */
// class SampleCostFunction2 : public GenericCostFunction
// {
// public:
//  /**
//   * We overload the operator as is required in GenericCostFunction to
//   * calculate the cost of the function.
//   *
//   * Behaviour:
//   *   cost = y^2 + 4xy +x^3 - 2x
//   *
//   * Intended use :
//   *   SampleCostFunction2 Func;
//   *   cost = Func(const std::vector<double> x);
//   *
//   * NOTE that the returned cost is the sqare root of the evaluation
//   *      due to the fact the Nelder Mead class will square the cost.
//   *
//   * Parameter:
//   *   \param x (input) vector containing coordinates of x and y
//   *   \return cost
//   */
//  std::vector<double> operator()(const std::vector<double> &x)
//  {
//    std::vector<double> ret {0};
//    ret[0] = sqrt(fabs(x[1]*x[1] + 4*x[1]*x[0] + x[0]*x[0]*x[0]- 2*x[0]));
//    return ret;
//  }
// };


class SampleCostFunction3 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = sqrt(sqrt(x[1]*x[1] + x[0]*x[0]));
    return ret;
  }
};

// /**
// * The Sample function 4 is defined as
// *
// *   (x^2 + (y - 11)^2 + (y^2 - 7 + x)^2).
// *
// * Number of dimensions = 2
// * The global minimum is: (3, 2)
// * Initial guess: (0, 0)
// *
// * Reference:
// */
// class SampleCostFunction4 : public GenericCostFunction
// {
// public:
//  /**
//   * We overload the operator as is required in GenericCostFunction to
//   * calculate the cost of the function.
//   *
//   * Behaviour:
//   *   cost = x^2 + (y - 11)^2 + (y^2 - 7 + x)^2
//   *
//   * Intended use :
//   *   SampleCostFunction4 Func;
//   *   cost = Func(const std::vector<double> x);
//   *
//   * NOTE that the returned cost is the sqare root of the evaluation
//   *      due to the fact the Nelder Mead class will square the cost.
//   *
//   * Parameters:
//   *   \param x (input) vector containing coordinates of x and y
//   *   \return cost
//   */
//  std::vector<double> operator()(const std::vector<double> &x)
//  {
//    std::vector<double> ret {0};
//    ret[0] = sqrt(fabs(x[0]*x[0] + x[1] - 11)*(x[0]*x[0] + x[1] - 11)
//      + (x[1]*x[1] - 7 + x[0])*(x[1]*x[1] - 7 + x[0]));
//    return ret;
//  }
// };

class SampleCostFunction5 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = sqrt(x[0]*x[0] + x[1]*x[1] + x[2]*x[2]);
    return ret;
  }
};

class SampleCostFunction6 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = sqrt(x[0]*x[0] + x[1]*x[1] + x[2]*x[2] + x[3]*x[3] + x[4]*x[4]);
    return ret;
  }
};

class SampleCostFunction7 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = sqrt(fabs(1/((x[0] - 1.05)*x[1]*x[2])));
    return ret;
  }
};

class SampleCostFunction8 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = sqrt(1/(x[0]*x[0]) + 1/(x[1]*x[1]));
    return ret;
  }
};

class SampleCostFunction9 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = sqrt(fabs(1/((x[0] - 0.95)*x[1]*x[2])));
    return ret;
  }
};

class Asymptote2D : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = (fabs(-1/(x[0]*x[1] + x[0]*(x[0] - 3))));
    return ret;
  }
};

class FailFindMin : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    if (x[0] < 1 && x[1]>1.025) {
      ret[0] = 1/(x[0]-x[0]);
    }
    else {
      ret[0] = 5*x[0]*x[0] + x[1]*x[1];
    }
    ret[0] = sqrt(ret[0]);
    return ret;
  }
};

class FailFindMin2 : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = (x[0]-90)*(x[0]-90) + (x[1]+50)*(x[1]+50) + 1/(x[0]-102);
    if (x[0] < 105) {
      ret[0] += 1/(1 - (0.5*x[0]+x[1]+48.5));
    }

    ret[0] = sqrt(fabs(ret[0]));
    return ret;
  }
};

class PowellSingular : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    std::vector<double> ret {0};
    ret[0] = ((x[0] + 10*x[1])*(x[0] + 10*x[1]) + 5*(x[2] -
     x[3])*(x[2] - x[3]) + (x[1] - 2*x[2])*(x[1] - 2*x[2])*(x[1] -
     2*x[2])*(x[1] - 2*x[2]) + 10*(x[0] - x[3])*(x[0] -
     x[3])*(x[0] - x[3])*(x[0] - x[3]));

    ret[0] = sqrt(fabs(ret[0]));
    return ret;
  }
};

class AlwaysInfinite : public GenericCostFunction
{
 public:
  std::vector<double> operator()(const std::vector<double> &x)
  {
    const double temp = (x[0] >= 1.0) ? x[0] : 1.0;
    std::vector<double> ret {0};
    ret[0] = (std::numeric_limits<double>::max()*1000.0*temp);
    return ret;
  }
};
}  // UnitTests namespace
}  // Unfit namespace
#endif