Golden section search
(include/search/golden_section_search.hpp)

• View this file on GitHub
• Last update: 2021-08-11 13:32:54+09:00
• Include: #include "include/search/golden_section_search.hpp"

黄金分割探索によって凸関数の区間内の最小値(または最大値)とそれを達成する引数を求める関数が定義されています。

---

golden_section_search<true>(a, b, f, diff)

区間 $I = [a, b]$ で定義される凸関数 $f: I \to \mathbf{R}$ に対して $\mathrm{argmin}{x \in I} \, f(x)$ と $\mathrm{min}{x \in I} f(x)$ をこの順番で持つ std::pair 型の値を返します。

diff は $\mathrm{argmin} \, f(x)$ に対する許容誤差です。

golden_section_search<false>(a, b, f, diff)

区間 $I = [a, b]$ で定義される凸関数 $f: I \to \mathbf{R}$ に対して $\mathrm{argmax}{x \in I} \, f(x)$ と $\mathrm{max}{x \in I} f(x)$ をこの順番で持つ std::pair 型の値を返します。

diff は $\mathrm{argmax} \, f(x)$ に対する許容誤差です。

---

Verified with

• test/search/golden_section_search/1.test.cpp
• test/search/golden_section_search/2.test.cpp

Code

//! @file golden_section_search.hpp

#ifndef CP_LIBRARY_GOLDEN_SECTION_SEARCH_HPP
#define CP_LIBRARY_GOLDEN_SECTION_SEARCH_HPP

#include <cassert>
#include <cmath>
#include <type_traits>
#include <utility>

#ifndef CP_LIBRARY_ASSERT
//! @brief Assert macro
#  define CP_LIBRARY_ASSERT(...) assert(__VA_ARGS__)
#  define CP_LIBRARY_ASSERT_NOT_DEFINED
#endif

namespace lib {

//! @brief Function to find the maximum or minimum value of a convex function f(x) when x is in [a, b]
//! @tparam minimize Set to true if you want to minimize the f(x), otherwise set to false.
//! @tparam RealType type of x (deduced from parameters)
//! @tparam Func type of f (deduced from parameters)
//! @param low lower bound (a)
//! @param high upper bound (b)
//! @param f function to minimize or maximize
//! @param diff acceptable error
//! @return std::pair { argmin(f(x)), min(f(x)) } (or argmax & max)
//! @note time complexity: O(log((high - low) / diff * (time needed to calculate f(x))))
template <bool minimize, typename RealType, typename Func>
[[nodiscard]] auto golden_section_search(RealType low, RealType high, const Func& f, const RealType diff = 1e-9L) {
  using F_ResType = decltype(f(std::declval<RealType>()));
  CP_LIBRARY_ASSERT(low <= high);

  using std::sqrt;
  const RealType phi        = (1 + sqrt(RealType(5.0L))) / 2;
  const RealType phi_plus_1 = phi + 1;

  RealType mid_low     = (low * phi + high) / phi_plus_1;
  RealType mid_high    = (low + high * phi) / phi_plus_1;
  F_ResType score_low  = f(mid_low);
  F_ResType score_high = f(mid_high);

  while (high - low > diff) {
    if (minimize ^ (score_low < score_high)) {
      low        = mid_low;
      mid_low    = mid_high;
      score_low  = score_high;
      mid_high   = (low + high * phi) / phi_plus_1;
      score_high = f(mid_high);
    } else {
      high       = mid_high;
      mid_high   = mid_low;
      score_high = score_low;
      mid_low    = (low * phi + high) / phi_plus_1;
      score_low  = f(mid_low);
    }
  }

  return std::pair {low, score_low};
}

}  // namespace lib

#ifdef CP_LIBRARY_ASSERT_NOT_DEFINED
#  undef CP_LIBRARY_ASSERT
#  undef CP_LIBRARY_ASSERT_NOT_DEFINED
#endif

#endif  // CP_LIBRARY_GOLDEN_SECTION_SEARCH_HPP

#line 1 "include/search/golden_section_search.hpp"

//! @file golden_section_search.hpp

#ifndef CP_LIBRARY_GOLDEN_SECTION_SEARCH_HPP
#define CP_LIBRARY_GOLDEN_SECTION_SEARCH_HPP

#include <cassert>
#include <cmath>
#include <type_traits>
#include <utility>

#ifndef CP_LIBRARY_ASSERT
//! @brief Assert macro
#  define CP_LIBRARY_ASSERT(...) assert(__VA_ARGS__)
#  define CP_LIBRARY_ASSERT_NOT_DEFINED
#endif

namespace lib {

//! @brief Function to find the maximum or minimum value of a convex function f(x) when x is in [a, b]
//! @tparam minimize Set to true if you want to minimize the f(x), otherwise set to false.
//! @tparam RealType type of x (deduced from parameters)
//! @tparam Func type of f (deduced from parameters)
//! @param low lower bound (a)
//! @param high upper bound (b)
//! @param f function to minimize or maximize
//! @param diff acceptable error
//! @return std::pair { argmin(f(x)), min(f(x)) } (or argmax & max)
//! @note time complexity: O(log((high - low) / diff * (time needed to calculate f(x))))
template <bool minimize, typename RealType, typename Func>
[[nodiscard]] auto golden_section_search(RealType low, RealType high, const Func& f, const RealType diff = 1e-9L) {
  using F_ResType = decltype(f(std::declval<RealType>()));
  CP_LIBRARY_ASSERT(low <= high);

  using std::sqrt;
  const RealType phi        = (1 + sqrt(RealType(5.0L))) / 2;
  const RealType phi_plus_1 = phi + 1;

  RealType mid_low     = (low * phi + high) / phi_plus_1;
  RealType mid_high    = (low + high * phi) / phi_plus_1;
  F_ResType score_low  = f(mid_low);
  F_ResType score_high = f(mid_high);

  while (high - low > diff) {
    if (minimize ^ (score_low < score_high)) {
      low        = mid_low;
      mid_low    = mid_high;
      score_low  = score_high;
      mid_high   = (low + high * phi) / phi_plus_1;
      score_high = f(mid_high);
    } else {
      high       = mid_high;
      mid_high   = mid_low;
      score_high = score_low;
      mid_low    = (low * phi + high) / phi_plus_1;
      score_low  = f(mid_low);
    }
  }

  return std::pair {low, score_low};
}

}  // namespace lib

#ifdef CP_LIBRARY_ASSERT_NOT_DEFINED
#  undef CP_LIBRARY_ASSERT
#  undef CP_LIBRARY_ASSERT_NOT_DEFINED
#endif

#endif  // CP_LIBRARY_GOLDEN_SECTION_SEARCH_HPP

Back to top page