// Copyright (C) 2017 Davis E. King (davis@dlib.net)
// License: Boost Software License See LICENSE.txt for the full license.
#include <dlib/global_optimization.h>
#include <dlib/statistics.h>
#include <sstream>
#include <string>
#include <cstdlib>
#include <ctime>
#include <vector>
#include <dlib/rand.h>
#include "tester.h"
namespace
{
using namespace test;
using namespace dlib;
using namespace std;
logger dlog("test.global_optimization");
// ----------------------------------------------------------------------------------------
void test_upper_bound_function(double relative_noise_magnitude, double solver_eps)
{
print_spinner();
dlog << LINFO << "test_upper_bound_function, relative_noise_magnitude="<< relative_noise_magnitude << ", solver_eps=" << solver_eps;
auto rosen = [](const matrix<double,0,1>& x) { return -1*( 100*std::pow(x(1) - x(0)*x(0),2.0) + std::pow(1 - x(0),2)); };
dlib::rand rnd;
auto make_rnd = [&rnd]() { matrix<double,0,1> x(2); x = 2*rnd.get_random_double(), 2*rnd.get_random_double(); return x; };
std::vector<function_evaluation> evals;
for (int i = 0; i < 100; ++i)
{
auto x = make_rnd();
evals.emplace_back(x,rosen(x));
}
upper_bound_function ub(evals, relative_noise_magnitude, solver_eps);
DLIB_TEST(ub.num_points() == (long)evals.size());
DLIB_TEST(ub.dimensionality() == 2);
for (auto& ev : evals)
{
dlog << LINFO << ub(ev.x) - ev.y;
DLIB_TEST_MSG(ub(ev.x) - ev.y > -1e10, ub(ev.x) - ev.y);
}
for (int i = 0; i < 100; ++i)
{
auto x = make_rnd();
evals.emplace_back(x,rosen(x));
ub.add(evals.back());
}
DLIB_TEST(ub.num_points() == (long)evals.size());
DLIB_TEST(ub.dimensionality() == 2);
for (auto& ev : evals)
{
dlog << LINFO << ub(ev.x) - ev.y;
DLIB_TEST_MSG(ub(ev.x) - ev.y > -1e10, ub(ev.x) - ev.y);
}
if (solver_eps < 0.001)
{
dlog << LINFO << "out of sample points: ";
for (int i = 0; i < 10; ++i)
{
auto x = make_rnd();
dlog << LINFO << ub(x) - rosen(x);
DLIB_TEST_MSG(ub(x) - rosen(x) > 1e-10, ub(x) - rosen(x));
}
}
}
// ----------------------------------------------------------------------------------------
double complex_holder_table ( double x0, double x1)
{
// The regular HolderTable function
//return -std::abs(sin(x0)*cos(x1)*exp(std::abs(1-std::sqrt(x0*x0+x1*x1)/pi)));
// My more complex version of it with discontinuities and more local minima.
double sign = 1;
for (double j = -4; j < 9; j += 0.5)
{
if (j < x0 && x0 < j+0.5)
x0 += sign*0.25;
sign *= -1;
}
// HolderTable function tilted towards 10,10
return -std::abs(sin(x0)*cos(x1)*exp(std::abs(1-std::sqrt(x0*x0+x1*x1)/pi))) +(x0+x1)/10 + sin(x0*10)*cos(x1*10);
}
// ----------------------------------------------------------------------------------------
void test_global_function_search()
{
function_spec spec{{-10,-10}, {10,10}};
function_spec spec2{{-10,-10, -50}, {10,10, 50}};
global_function_search opt({spec, spec, spec2});
dlib::rand rnd;
bool found_optimal_point = false;
for (int i = 0; i < 400 && !found_optimal_point; ++i)
{
print_spinner();
std::vector<function_evaluation_request> nexts;
for (int k = 0; k < rnd.get_integer_in_range(1,4); ++k)
nexts.emplace_back(opt.get_next_x());
for (auto& next : nexts)
{
switch (next.function_idx())
{
case 0: next.set( -complex_holder_table(next.x()(0), next.x()(1))); break;
case 1: next.set( -10*complex_holder_table(next.x()(0), next.x()(1))); break;
case 2: next.set( -2*complex_holder_table(next.x()(0), next.x()(1))); break;
default: DLIB_TEST(false); break;
}
matrix<double,0,1> x;
double y;
size_t function_idx;
opt.get_best_function_eval(x,y,function_idx);
/*
cout << "\ni: "<< i << endl;
cout << "best eval x: "<< trans(x);
cout << "best eval y: "<< y << endl;
cout << "best eval function index: "<< function_idx << endl;
*/
if (std::abs(y - 10*21.9210397) < 0.0001)
{
found_optimal_point = true;
break;
}
}
}
DLIB_TEST(found_optimal_point);
}
// ----------------------------------------------------------------------------------------
void test_find_max_global(
)
{
print_spinner();
auto rosen = [](const matrix<double,0,1>& x) { return -1*( 100*std::pow(x(1) - x(0)*x(0),2.0) + std::pow(1 - x(0),2)); };
auto result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, max_function_calls(100), 0);
matrix<double,0,1> true_x = {1,1};
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, max_function_calls(100));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, max_function_calls(100), 0, std::vector<function_evaluation>{function_evaluation({1.1, 0.9}, rosen({1.1, 0.9}))});
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, std::chrono::seconds(5));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, std::chrono::seconds(5), 0.0);
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, std::chrono::seconds(5), max_function_calls(100));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, max_function_calls(100), std::chrono::seconds(5));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, {false,false}, max_function_calls(100));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {2, 2}, {false,false}, max_function_calls(100), 0.0);
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global(rosen, {0.1, 0.1}, {0.9, 0.9}, {false,false}, max_function_calls(140));
true_x = {0.9, 0.81};
dlog << LINFO << "rosen, bounded at 0.9: " << trans(result.x);
DLIB_TEST_MSG(max(abs(true_x-result.x)) < 1e-5, max(abs(true_x-result.x)));
print_spinner();
result = find_max_global([](double x){ return -std::pow(x-2,2.0); }, -10, 10, max_function_calls(10), 0);
dlog << LINFO << "(x-2)^2: " << trans(result.x);
DLIB_TEST(result.x.size()==1);
DLIB_TEST(std::abs(result.x - 2) < 1e-9);
print_spinner();
unsigned int normal_evals=0, early_evals=0;
auto normal_result = find_max_global([&normal_evals](double x){ normal_evals++; return -std::pow(x-2,2.0); }, -10, 1, max_function_calls(10));
dlog << LINFO << "(x-2)^2, bound at 1: " << trans(normal_result.x);
DLIB_TEST(normal_result.x.size()==1);
DLIB_TEST(std::abs(normal_result.x - 1) < 1e-9);
print_spinner();
result = find_max_global([](double x){ return -std::pow(x-2,2.0); }, -10, 1, std::chrono::seconds(2));
dlog << LINFO << "(x-2)^2, bound at 1: " << trans(result.x);
DLIB_TEST(result.x.size()==1);
DLIB_TEST(std::abs(result.x - 1) < 1e-9);
print_spinner();
constexpr auto close_enough = -16.0;
auto early_result = find_max_global([&early_evals](double x){ early_evals++; return -std::pow(x-2,2.0); }, -10, 1, max_function_calls(10), 0.0, std::vector<function_evaluation>{}, [&](double y){ return (y >= close_enough);});
dlog << LINFO << "(x-2)^2, bound at 1: " << trans(early_result.x);
DLIB_TEST(early_result.x.size()==1);
DLIB_TEST(std::abs(early_result.y) <= std::abs(close_enough));
DLIB_TEST(std::abs(early_result.x - 1) <= 4);
DLIB_TEST(normal_evals >= early_evals);
print_spinner();
result = find_max_global([](double a, double b){ return -complex_holder_table(a,b);},
{-10, -10}, {10, 10}, max_function_calls(400), 0);
dlog << LINFO << "complex_holder_table y: "<< result.y;
DLIB_TEST_MSG(std::abs(result.y - 21.9210397) < 0.0001, std::abs(result.y - 21.9210397));
}
// ----------------------------------------------------------------------------------------
void test_find_min_global(
)
{
print_spinner();
auto rosen = [](const matrix<double,0,1>& x) { return +1*( 100*std::pow(x(1) - x(0)*x(0),2.0) + std::pow(1 - x(0),2)); };
auto result = find_min_global(rosen, {0.1, 0.1}, {2, 2}, max_function_calls(100), 0);
matrix<double,0,1> true_x = {1,1};
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(min(abs(true_x-result.x)) < 1e-5, min(abs(true_x-result.x)));
print_spinner();
result = find_min_global(rosen, {0.1, 0.1}, {2, 2}, max_function_calls(100));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(min(abs(true_x-result.x)) < 1e-5, min(abs(true_x-result.x)));
print_spinner();
result = find_min_global(rosen, {0.1, 0.1}, {2, 2}, std::chrono::seconds(5));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(min(abs(true_x-result.x)) < 1e-5, min(abs(true_x-result.x)));
print_spinner();
result = find_min_global(rosen, {0.1, 0.1}, {2, 2}, {false,false}, max_function_calls(100));
dlog << LINFO << "rosen: " << trans(result.x);
DLIB_TEST_MSG(min(abs(true_x-result.x)) < 1e-5, min(abs(true_x-result.x)));
print_spinner();
result = find_min_global(rosen, {0, 0}, {0.9, 0.9}, {false,false}, max_function_calls(100));
true_x = {0.9, 0.81};
dlog << LINFO << "rosen, bounded at 0.9: " << trans(result.x);
DLIB_TEST_MSG(min(abs(true_x-result.x)) < 1e-5, min(abs(true_x-result.x)));
print_spinner();
result = find_min_global([](double x){ return std::pow(x-2,2.0); }, -10, 10, max_function_calls(10), 0);
dlog << LINFO << "(x-2)^2: " << trans(result.x);
DLIB_TEST(result.x.size()==1);
DLIB_TEST(std::abs(result.x - 2) < 1e-9);
print_spinner();
result = find_min_global([](double x){ return std::pow(x-2,2.0); }, -10, 1, max_function_calls(10));
dlog << LINFO << "(x-2)^2, bound at 1: " << trans(result.x);
DLIB_TEST(result.x.size()==1);
DLIB_TEST(std::abs(result.x - 1) < 1e-9);
print_spinner();
result = find_min_global([](double x){ return std::pow(x-2,2.0); }, -10, 1, std::chrono::seconds(2));
dlog << LINFO << "(x-2)^2, bound at 1: " << trans(result.x);
DLIB_TEST(result.x.size()==1);
DLIB_TEST(std::abs(result.x - 1) < 1e-9);
print_spinner();
result = find_min_global([](double a, double b){ return complex_holder_table(a,b);},
{-10, -10}, {10, 10}, max_function_calls(400), 0);
dlog << LINFO << "complex_holder_table y: "<< result.y;
DLIB_TEST_MSG(std::abs(result.y + 21.9210397) < 0.0001, std::abs(result.y + 21.9210397));
}
// ----------------------------------------------------------------------------------------
class global_optimization_tester : public tester
{
public:
global_optimization_tester (
) :
tester ("test_global_optimization",
"Runs tests on the global optimization components.")
{}
void perform_test (
)
{
test_upper_bound_function(0.01, 1e-6);
test_upper_bound_function(0.0, 1e-6);
test_upper_bound_function(0.0, 1e-1);
test_global_function_search();
test_find_max_global();
test_find_min_global();
}
} a;
}