diff --git a/src/problem/base_meta.h b/src/problem/base_meta.h index 7b4275df..1e68e1c9 100644 --- a/src/problem/base_meta.h +++ b/src/problem/base_meta.h @@ -68,6 +68,11 @@ class __PAGMO_VISIBLE base_meta : public base bool compare_fc_impl(const fitness_vector &f1, const constraint_vector &c1, const fitness_vector &f2, const constraint_vector &c2) const {return m_original_problem->compare_fc_impl(f1,c1,f2,c2);} private: + // Computing new optima from the optima of the original problem. + // Applicable to only some transformations. + virtual std::vector transform_x(const std::vector &) const + { return std::vector(); } + friend class boost::serialization::access; template void serialize(Archive &ar, const unsigned int) diff --git a/src/problem/himmelblau.cpp b/src/problem/himmelblau.cpp index f4934600..0017dca4 100644 --- a/src/problem/himmelblau.cpp +++ b/src/problem/himmelblau.cpp @@ -30,7 +30,11 @@ namespace pagmo { namespace problem { /// Default constructor. -himmelblau::himmelblau():base(-6.,6.,2) {} +himmelblau::himmelblau():base(-6.,6.,2) +{ + // initialize best solution + initialize_best(); +} /// Clone method. base_ptr himmelblau::clone() const @@ -51,6 +55,30 @@ std::string himmelblau::get_name() const return "Himmelblau"; } +/// Initialize the Himmelblau's function four local minima. +void himmelblau::initialize_best(void) +{ + const int min_count = 4; + const int x_dimension = 2; + + const double x_vector[][x_dimension] = { + {3.0, 2.0}, + {-2.805118, 3.131312}, + {-3.779310, -3.283186}, + {3.584428, -1.848126} + }; + + std::vector best_x(min_count); + + for (int i = 0; i < min_count; ++i) { + decision_vector x(x_dimension); + std::copy(x_vector[i],x_vector[i] + x_dimension,x.begin()); + best_x[i] = x; + } + + set_best_x(best_x); +} + }} BOOST_CLASS_EXPORT_IMPLEMENT(pagmo::problem::himmelblau) diff --git a/src/problem/himmelblau.h b/src/problem/himmelblau.h index 2c628d6c..d550e10d 100644 --- a/src/problem/himmelblau.h +++ b/src/problem/himmelblau.h @@ -61,6 +61,8 @@ class __PAGMO_VISIBLE himmelblau: public base protected: void objfun_impl(fitness_vector &, const decision_vector &) const; private: + void initialize_best(void); + friend class boost::serialization::access; template void serialize(Archive &ar, const unsigned int) diff --git a/src/problem/normalized.cpp b/src/problem/normalized.cpp index 8afec8f5..13bf725d 100644 --- a/src/problem/normalized.cpp +++ b/src/problem/normalized.cpp @@ -53,6 +53,8 @@ normalized::normalized(const base & p): m_normalization_scale(p.get_dimension(),0) { configure_new_bounds(); + std::vector new_best_x = transform_x(p.get_best_x()); + set_best_x(new_best_x); } /// Clone method. @@ -73,6 +75,25 @@ void normalized::configure_new_bounds() set_bounds(-1, 1); } +/// Compute normalized vectors for the meta-problem from the original problem. +/* + * @param[in] x vectors of the original problem + * @param[out] vectors x normalized + */ +std::vector normalized::transform_x(const std::vector &x) const +{ + const base::size_type cnt = x.size(); + std::vector new_x = x; + + for (base::size_type i = 0; i < cnt; ++i) { + for (base::size_type j = 0; j < x[i].size(); ++j) { + new_x[i][j] = (x[i][j] - m_normalization_center[j]) / m_normalization_scale[j]; + } + } + + return new_x; +} + /// Returns the de-normalized version of the decision variables decision_vector normalized::denormalize(const decision_vector& x) const { diff --git a/src/problem/normalized.h b/src/problem/normalized.h index d30cb0ff..decf6d73 100644 --- a/src/problem/normalized.h +++ b/src/problem/normalized.h @@ -58,6 +58,7 @@ class __PAGMO_VISIBLE normalized : public base_meta void compute_constraints_impl(constraint_vector &, const decision_vector &) const; private: void configure_new_bounds(); + std::vector transform_x(const std::vector &) const; friend class boost::serialization::access; template diff --git a/src/problem/rotated.cpp b/src/problem/rotated.cpp index 54803880..812a20b6 100644 --- a/src/problem/rotated.cpp +++ b/src/problem/rotated.cpp @@ -63,6 +63,8 @@ rotated::rotated(const base &p, const Eigen::MatrixXd &rotation ): pagmo_throw(value_error,"Input problem has an integer dimension. Cannot rotate it."); } configure_new_bounds(); + std::vector new_best_x = transform_x(p.get_best_x()); + set_best_x(new_best_x); } /** @@ -107,6 +109,8 @@ rotated::rotated(const base &p, pagmo_throw(value_error,"The input matrix seems not to be orthonormal (to a tolerance of 1e-5)"); } configure_new_bounds(); + std::vector new_best_x = transform_x(p.get_best_x()); + set_best_x(new_best_x); } /** @@ -138,6 +142,8 @@ rotated::rotated(const base &p): pagmo_throw(value_error,"Input problem has an integer dimension. Cannot rotate it."); } configure_new_bounds(); + std::vector new_best_x = transform_x(p.get_best_x()); + set_best_x(new_best_x); } /// Clone method. @@ -170,6 +176,39 @@ void rotated::configure_new_bounds() set_bounds(-sqrt(2), sqrt(2)); } +/// Compute rotated vectors for the meta-problem from the original problem. +/* + * @param[in] x vectors of the original problem + * @param[out] vectors x rotated + */ +std::vector rotated::transform_x(const std::vector &x) const +{ + const base::size_type cnt = x.size(); + std::vector new_x = x; + + for (base::size_type i = 0; i < cnt; ++i) { + // 1. normalize the vector x[i] + decision_vector x_normed(x[i].size()); + for (base::size_type j = 0; j < x[i].size(); ++j) { + x_normed[j] = (x[i][j] - m_normalize_translation[j]) / m_normalize_scale[j]; + } + + // 2. rotate the normalized vector + Eigen::VectorXd x_normed_vec = Eigen::VectorXd::Zero(x_normed.size()); + for (base::size_type j = 0; j < x_normed.size(); ++j) { + x_normed_vec(j) = x_normed[j]; + } + Eigen::VectorXd x_rotated_vec = m_Rotate * x_normed_vec; + + // Store the normalized and rotated vector + for (base::size_type j = 0; j < x[i].size(); ++j) { + new_x[i][j] = x_rotated_vec(j); + } + } + + return new_x; +} + // Used to normalize the original upper and lower bounds // to [-1, 1], at each dimension decision_vector rotated::normalize_to_center(const decision_vector& x) const diff --git a/src/problem/rotated.h b/src/problem/rotated.h index ff3120ba..3b1301a5 100644 --- a/src/problem/rotated.h +++ b/src/problem/rotated.h @@ -64,6 +64,7 @@ class __PAGMO_VISIBLE rotated : public base_meta private: void configure_new_bounds(); + std::vector transform_x(const std::vector &) const; decision_vector normalize_to_center(const decision_vector& x) const; decision_vector denormalize_to_original(const decision_vector& x) const; diff --git a/src/problem/scaled.cpp b/src/problem/scaled.cpp index 0cafe371..0b2ca6a4 100644 --- a/src/problem/scaled.cpp +++ b/src/problem/scaled.cpp @@ -102,6 +102,16 @@ void scaled::compute_constraints_impl(constraint_vector &c, const decision_vecto m_original_problem->compute_constraints(c, x); } +/// Compute transformed vectors for the meta-problem from the original problem with the fitness scaled. +/* + * @param[in] x vectors of the original problem + * @param[out] vectors of the scaled problem + */ +std::vector scaled::transform_x(const std::vector &x) const +{ + // Only the fitness is scaled, not the space itself. + return x; +} std::string scaled::get_name() const { diff --git a/src/problem/scaled.h b/src/problem/scaled.h index 2a1e92c7..d3087875 100644 --- a/src/problem/scaled.h +++ b/src/problem/scaled.h @@ -58,6 +58,7 @@ class __PAGMO_VISIBLE scaled : public base_meta void objfun_impl(fitness_vector &, const decision_vector &) const; void compute_constraints_impl(constraint_vector &, const decision_vector &) const; private: + std::vector transform_x(const std::vector &) const; friend class boost::serialization::access; template diff --git a/src/problem/shifted.cpp b/src/problem/shifted.cpp index 94f0dc70..10579041 100644 --- a/src/problem/shifted.cpp +++ b/src/problem/shifted.cpp @@ -57,6 +57,8 @@ shifted::shifted(const base & p, pagmo_throw(value_error,"The size of the shifting vector must be equal to the problem dimension"); } configure_shifted_bounds(m_translation); + std::vector new_best_x = transform_x(p.get_best_x()); + set_best_x(new_best_x); } @@ -82,6 +84,8 @@ shifted::shifted(const base & p, m_translation(decision_vector(p.get_dimension(), t)) { configure_shifted_bounds(m_translation); + std::vector new_best_x = transform_x(p.get_best_x()); + set_best_x(new_best_x); } /** @@ -107,6 +111,8 @@ shifted::shifted(const base & p): m_translation[i] = (2*((double) rand() / (RAND_MAX))-1) * (p.get_ub()[i]-p.get_lb()[i]); } configure_shifted_bounds(m_translation); + std::vector new_best_x = transform_x(p.get_best_x()); + set_best_x(new_best_x); } /// Clone method. @@ -159,6 +165,25 @@ void shifted::compute_constraints_impl(constraint_vector &c, const decision_vect m_original_problem->compute_constraints(c, x_translated); } +/// Compute shifted vectors for the meta-problem from the original problem. +/* + * @param[in] x vectors of the original problem + * @param[out] vectors x shifted + */ +std::vector shifted::transform_x(const std::vector &x) const +{ + const base::size_type cnt = x.size(); + std::vector new_x = x; + + for (base::size_type i = 0; i < cnt; ++i) { + for (base::size_type j = 0; j < x[i].size(); ++j) { + new_x[i][j] = x[i][j] + m_translation[j]; + } + } + + return new_x; +} + /** * Gets the shift vector which defines the problem * diff --git a/src/problem/shifted.h b/src/problem/shifted.h index 22096d1d..c43d187a 100644 --- a/src/problem/shifted.h +++ b/src/problem/shifted.h @@ -62,6 +62,7 @@ class __PAGMO_VISIBLE shifted : public base_meta void compute_constraints_impl(constraint_vector &, const decision_vector &) const; private: void configure_shifted_bounds(const decision_vector &); + std::vector transform_x(const std::vector &) const; friend class boost::serialization::access; template diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 0ac9bdb7..a8e55f3a 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -50,6 +50,10 @@ ADD_EXECUTABLE(test_tsp test_tsp.cpp) TARGET_LINK_LIBRARIES(test_tsp ${MANDATORY_LIBRARIES} pagmo_static) ADD_TEST(test_tsp test_tsp) +ADD_EXECUTABLE(test_optima_transformations test_optima_transformations.cpp) +TARGET_LINK_LIBRARIES(test_optima_transformations ${MANDATORY_LIBRARIES} pagmo_static) +ADD_TEST(test_optima_transformations test_optima_transformations) + ADD_EXECUTABLE(test_archipelago test_archipelago.cpp) TARGET_LINK_LIBRARIES(test_archipelago ${MANDATORY_LIBRARIES} pagmo_static) ADD_TEST(test_archipelago test_archipelago) diff --git a/tests/best_solutions_test.cpp b/tests/best_solutions_test.cpp index 42dbde79..7faaf253 100644 --- a/tests/best_solutions_test.cpp +++ b/tests/best_solutions_test.cpp @@ -180,7 +180,7 @@ int main() std::cout << " fitness passes, "; } else{ - std::cout << " fitness failed!"<get_best_c().at(j); if(is_eq(best_c_computed, best_c, EPS)){ - std::cout << " constraints passes."; + std::cout << " constraints passes." << std::endl; } else{ - std::cout << " constraints failed!"< +#include +#include +#include "../src/pagmo.h" + +using namespace pagmo; + +const double EPS = 10e-9; + +bool is_eq(const decision_vector & x1, const decision_vector & x2){ + if (x1.size() != x2.size()) return false; + for (unsigned int i = 0; i < x1.size(); ++i){ + if (fabs(x1[i] - x2[i]) > EPS) return false; + } + return true; +} + +// Going to test the code on the Himmelblau function, +// see http://en.wikipedia.org/wiki/Himmelblau%27s_function +// The four minima of the function in the range [-6,6]x[-6,6] are: +// (x0,y0) = (3.0, 2.0) +// (x1,y1) = (-2.805118,3.131312) +// (x2,y2) = (-3.779310,-3.283186) +// (x3,y3) = (3.584428,-1.848126) + +// Test Himmelblau function minima in the bounding box [-6,6]: +const std::vector best_x = { + {3.0, 2.0}, + {-2.805118,3.131312}, + {-3.779310,-3.283186}, + {3.584428,-1.848126}, +}; + +int test_shifted() +{ + problem::himmelblau prob; + problem::shifted nonzero_shifted_prob(prob, {100,50}); + // retrieve the shifted minima + const std::vector shifted_minima = nonzero_shifted_prob.get_best_x(); + + // just check that we have actualy moved (somewhere) + bool all_equal = true; + for (unsigned int i = 0; i < best_x.size(); ++i) { + if (!is_eq(best_x[i],shifted_minima[i])) all_equal = false; + } + if (all_equal) return 1; + + + // Perform a number of shifts that walk around and get back to the origin in the end. + // It is possible to add any combination of translations satisfying such a condition. + const std::vector> shifts = { + {1., 1.}, + {-2., 0.}, + {1., -1.}, + + {1., 1.}, + {-2., 0.}, + {1., -1.}, + }; + + // test the zero shift + problem::shifted zero_shifted_prob(prob, {0,0}); + const std::vector zero_shifted_minima = zero_shifted_prob.get_best_x(); + for (unsigned int i = 0; i < best_x.size(); ++i) { + if (!is_eq(best_x[i],zero_shifted_minima[i])) return 1; + } + + std::vector problems; + problems.push_back(zero_shifted_prob); + + for (unsigned int i = 0; i < shifts.size(); ++i) + problems.push_back(problem::shifted(problems.back(), shifts[i])); + + std::vector> expected_shifted_minima; + expected_shifted_minima.push_back(zero_shifted_prob.get_best_x()); + + for (unsigned int i = 0; i < shifts.size(); ++i) { + expected_shifted_minima.push_back({ + {expected_shifted_minima.back()[0][0] + shifts[i][0], expected_shifted_minima.back()[0][1] + shifts[i][1]}, + {expected_shifted_minima.back()[1][0] + shifts[i][0], expected_shifted_minima.back()[1][1] + shifts[i][1]}, + {expected_shifted_minima.back()[2][0] + shifts[i][0], expected_shifted_minima.back()[2][1] + shifts[i][1]}, + {expected_shifted_minima.back()[3][0] + shifts[i][0], expected_shifted_minima.back()[3][1] + shifts[i][1]}, + }); + } + + for (unsigned int i = 0; i < problems.size(); ++i) { + const std::vector shifted_best_x = problems[i].get_best_x(); + for (unsigned int j = 0; j < shifted_best_x.size(); ++j) { + if (!is_eq(expected_shifted_minima[i][j],shifted_best_x[j])) return 1; + } + } + + // Now check that we have actually returned to the origin back again. + for (unsigned int i = 0; i < best_x.size(); ++i) { + if (!is_eq(best_x[i],expected_shifted_minima.back()[i])) return 1; + } + + return 0; +} + +int test_rotated() +{ + // tak the original Himmelblau function + problem::himmelblau prob; + // and normalize it + problem::normalized normalized_prob(prob); + + const std::vector> rot_id_mat = {{1.,0.},{0.,1.}}; // identity matrix + problem::rotated rotated_id_prob(normalized_prob,rot_id_mat); + // Normalized Himmelblau rotated by the identity matrix and the corresponding minima are: + const std::vector rot_id_best_x = rotated_id_prob.get_best_x(); + + const std::vector expected_normalized_minima = { + {best_x[0][0] / 6, best_x[0][1] / 6}, + {best_x[1][0] / 6, best_x[1][1] / 6}, + {best_x[2][0] / 6, best_x[2][1] / 6}, + {best_x[3][0] / 6, best_x[3][1] / 6}, + }; + + for (unsigned int i = 0; i < expected_normalized_minima.size(); ++i) { + if (!is_eq(expected_normalized_minima[i], rot_id_best_x[i])) return 1; + } + + const double pi = std::acos(-1); + const double phi = pi / 2; + const std::vector> rot_mat90 = {{std::cos(phi), -std::sin(phi)},{std::sin(phi), std::cos(phi)}}; + + // Take the Himmelblau function and rotate by 90 degrees (plus normalize) + problem::rotated rotated_prob1(prob,rot_mat90); + + // Rotate again by 90 degrees (plus normalize) + problem::rotated rotated_prob2(rotated_prob1,rot_mat90); + + // Rotate again by 90 degrees (plus normalize) + problem::rotated rotated_prob3(rotated_prob2,rot_mat90); + + // Rotate again by 90 degrees (plus normalize). We get back to the original problem (normalized_prob), + // but there is a twist here, because along the way the bounds were normalized several + // times (exactly four time). Therefore, all the vectors are scaled by the factor (1 / sqrt(2)), + // except for the frist translation which is normalized by the factor 1 / 6. + // See src/problem/rotated.cpp for more details (function configure_new_bounds() in particular). + problem::rotated rotated_prob4(rotated_prob3,rot_mat90); + // Previous function rotated by 90 degrees counterclockwise and the minima are (after another normalization) + // basically the same as in the original problem, but normalized by the factor of (1 / sqrt(2)) ** 3: + const std::vector rot4_best_x = rotated_prob4.get_best_x(); + + const double sqrt2 = std::sqrt(2); + // First normalization goes from the bounding box of [-6,6]x[-6,6] to [-1,1]x[-1,1] + // and every other normalization goes from [-sqrt(2),sqrt(2)]x[-sqrt(2),sqrt(2)] to [-1,1]x[-1,1]. + const double coef = (1 / 6.) * (1 / sqrt2) * (1 / sqrt2) * (1 / sqrt2); + + const std::vector expected_shrunk_minima = { + {best_x[0][0] * coef, best_x[0][1] * coef}, + {best_x[1][0] * coef, best_x[1][1] * coef}, + {best_x[2][0] * coef, best_x[2][1] * coef}, + {best_x[3][0] * coef, best_x[3][1] * coef}, + }; + + for (unsigned int i = 0; i < expected_shrunk_minima.size(); ++i) { + if (!is_eq(expected_shrunk_minima[i], rot4_best_x[i])) return 1; + } + + return 0; +} + +int test_normalized() +{ + // tak the original Himmelblau function + problem::himmelblau prob; + // and normalize it + problem::normalized normalized_prob(prob); + // Normalizing bounds (to [-1,1]x[-1,1]) of the Himmelblau function and the corresponding transformed minima are: + const std::vector norm_best_x = normalized_prob.get_best_x(); + + const std::vector expected_normalized_minima = { + {best_x[0][0] / 6, best_x[0][1] / 6}, + {best_x[1][0] / 6, best_x[1][1] / 6}, + {best_x[2][0] / 6, best_x[2][1] / 6}, + {best_x[3][0] / 6, best_x[3][1] / 6}, + }; + + for (unsigned int i = 0; i < norm_best_x.size(); ++i) { + if (!is_eq(expected_normalized_minima[i], norm_best_x[i])) return 1; + } + + return 0; +} + +int main() +{ + return test_shifted() || test_rotated() || test_normalized(); +}