Skip to content

Commit

Permalink
A few changes
Browse files Browse the repository at this point in the history
  • Loading branch information
@mzhurovich
mzhurovich committed Mar 8, 2015
1 parent db559aa commit 3dac18f
Show file tree
Hide file tree
Showing 3 changed files with 43 additions and 38 deletions.
25 changes: 14 additions & 11 deletions fncas/fncas_mathutil.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,27 +31,30 @@

namespace fncas {

bool IsNormal(double arg) { return (std::isnormal(arg) || arg == 0.0); }
inline bool IsNormal(double arg) { return (std::isnormal(arg) || arg == 0.0); }

template <typename T>
inline typename std::enable_if<std::is_arithmetic<T>::value, T>::type DotProduct(const std::vector<T> &v1,
const std::vector<T> &v2) {
inline typename std::enable_if<std::is_arithmetic<T>::value, T>::type DotProduct(const std::vector<T>& v1,
const std::vector<T>& v2) {
#ifndef NDEBUG
assert(v1.size() == v2.size());
#endif
return std::inner_product(std::begin(v1), std::end(v1), std::begin(v2), static_cast<T>(0));
}

template <typename T>
inline typename std::enable_if<std::is_arithmetic<T>::value, T>::type L2Norm(const std::vector<T> &v) {
inline typename std::enable_if<std::is_arithmetic<T>::value, T>::type L2Norm(const std::vector<T>& v) {
return DotProduct(v, v);
}

template <typename T, typename = typename std::enable_if<std::is_arithmetic<T>::value, T>::type>
inline void FlipSign(std::vector<T> &v) {
inline void FlipSign(std::vector<T>& v) {
std::transform(std::begin(v), std::end(v), std::begin(v), std::negate<T>());
}

// Polak-Ribiere formula for conjugate gradient method
// http://en.wikipedia.org/wiki/Nonlinear_conjugate_gradient_method
inline double PolakRibiere(const std::vector<double> &g, const std::vector<double> &g_prev) {
inline double PolakRibiere(const std::vector<double>& g, const std::vector<double>& g_prev) {
const double beta = (L2Norm(g) - DotProduct(g, g_prev)) / L2Norm(g_prev);
if (IsNormal(beta)) {
return beta;
Expand All @@ -63,12 +66,12 @@ inline double PolakRibiere(const std::vector<double> &g, const std::vector<doubl
// Simplified backtracking line search algorithm with limited number of steps.
// Starts in `current_point` and searches for minimum in `direction`
// sequentially shrinking the step size. Returns new optimal point.
// Algorithm parameters: 0 < alpha < 1, 0 < beta < 1
// Algorithm parameters: 0 < alpha < 1, 0 < beta < 1.
template <class F, class G>
inline std::vector<double> BackTracking(F &&eval_function,
G &&eval_gradient,
const std::vector<double> &current_point,
const std::vector<double> &direction,
inline std::vector<double> BackTracking(F&& eval_function,
G&& eval_gradient,
const std::vector<double>& current_point,
const std::vector<double>& direction,
const double alpha = 0.5,
const double beta = 0.8,
const size_t max_steps = 100) {
Expand Down
40 changes: 21 additions & 19 deletions fncas/fncas_optimize.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -67,35 +67,38 @@ const T OptimizerParameters::GetValue(std::string name, T default_value) {
}

// TODO(mzhurovich+dkorolev): Make different implementation to use functors.
// Naive gradient descent that tries 3 different step sizes in each iteration.
// Searches for a local minimum of `F::Compute` function.
template <class F>
class GradientDescentOptimizer : noncopyable {
public:
GradientDescentOptimizer() {}
GradientDescentOptimizer(OptimizerParameters &params) {
GradientDescentOptimizer(OptimizerParameters& params) {
max_steps_ = params.GetValue("max_steps", max_steps_);
step_factor_ = params.GetValue("step_factor", step_factor_);
}
OptimizationResult Optimize(const std::vector<double> &starting_point);
OptimizationResult Optimize(const std::vector<double>& starting_point);

private:
size_t max_steps_ = 5000; // Maximum number of optimization steps.
double step_factor_ = 1.0; // Gradient is multiplied by this factor.
};

// Simple gradient descent optimizer using backtracking algorithm.
// Simple gradient descent optimizer with backtracking algorithm.
// Searches for a local minimum of `F::Compute` function.
template <class F>
class GradientDescentOptimizerBT : noncopyable {
public:
GradientDescentOptimizerBT() {}
GradientDescentOptimizerBT(OptimizerParameters &params) {
GradientDescentOptimizerBT(OptimizerParameters& params) {
min_steps_ = params.GetValue("min_steps", min_steps_);
max_steps_ = params.GetValue("max_steps", max_steps_);
bt_alpha_ = params.GetValue("bt_alpha", bt_alpha_);
bt_beta_ = params.GetValue("bt_beta", bt_beta_);
bt_max_steps_ = params.GetValue("bt_max_steps", bt_max_steps_);
grad_eps_ = params.GetValue("grad_eps", grad_eps_);
}
OptimizationResult Optimize(const std::vector<double> &starting_point);
OptimizationResult Optimize(const std::vector<double>& starting_point);

private:
size_t min_steps_ = 3; // Minimum number of optimization steps (ignoring early stopping).
Expand All @@ -107,20 +110,20 @@ class GradientDescentOptimizerBT : noncopyable {
};

// Optimizer that uses a combination of conjugate gradient method and
// backtracking line search.
// backtracking line search to find a local minimum of `F::Compute` function.
template <class F>
class ConjugateGradientOptimizer : noncopyable {
public:
ConjugateGradientOptimizer() {}
ConjugateGradientOptimizer(OptimizerParameters &params) {
ConjugateGradientOptimizer(OptimizerParameters& params) {
min_steps_ = params.GetValue("min_steps", min_steps_);
max_steps_ = params.GetValue("max_steps", max_steps_);
bt_alpha_ = params.GetValue("bt_alpha", bt_alpha_);
bt_beta_ = params.GetValue("bt_beta", bt_beta_);
bt_max_steps_ = params.GetValue("bt_max_steps", bt_max_steps_);
grad_eps_ = params.GetValue("grad_eps", grad_eps_);
}
OptimizationResult Optimize(const std::vector<double> &starting_point);
OptimizationResult Optimize(const std::vector<double>& starting_point);

private:
size_t min_steps_ = 3; // Minimum number of optimization steps (ignoring early stopping).
Expand All @@ -132,7 +135,7 @@ class ConjugateGradientOptimizer : noncopyable {
};

template <class F>
OptimizationResult GradientDescentOptimizer<F>::Optimize(const std::vector<double> &starting_point) {
OptimizationResult GradientDescentOptimizer<F>::Optimize(const std::vector<double>& starting_point) {
fncas::reset_internals_singleton();
const size_t dim = starting_point.size();
const fncas::x gradient_helper(dim);
Expand All @@ -159,7 +162,7 @@ OptimizationResult GradientDescentOptimizer<F>::Optimize(const std::vector<doubl
}

template <class F>
OptimizationResult GradientDescentOptimizerBT<F>::Optimize(const std::vector<double> &starting_point) {
OptimizationResult GradientDescentOptimizerBT<F>::Optimize(const std::vector<double>& starting_point) {
fncas::reset_internals_singleton();
const size_t dim = starting_point.size();
const fncas::x gradient_helper(dim);
Expand All @@ -169,11 +172,10 @@ OptimizationResult GradientDescentOptimizerBT<F>::Optimize(const std::vector<dou

for (size_t iteration = 0; iteration < max_steps_; ++iteration) {
auto direction = gi(current_point).gradient;
fncas::FlipSign(direction); // Direction = negative gradient.
auto new_point = BackTracking(fi, gi, current_point, direction, bt_alpha_, bt_beta_, bt_max_steps_);
current_point = new_point;
fncas::FlipSign(direction); // Going against the gradient to minimize the function.
current_point = BackTracking(fi, gi, current_point, direction, bt_alpha_, bt_beta_, bt_max_steps_);

// Simple early stopping by the norm of gradient.
// Simple early stopping by the norm of the gradient.
if (std::sqrt(fncas::L2Norm(direction)) < grad_eps_ && iteration >= min_steps_) {
break;
}
Expand All @@ -184,7 +186,7 @@ OptimizationResult GradientDescentOptimizerBT<F>::Optimize(const std::vector<dou

// TODO(mzhurovich): Implement more sophisticated version.
template <class F>
OptimizationResult ConjugateGradientOptimizer<F>::Optimize(const std::vector<double> &starting_point) {
OptimizationResult ConjugateGradientOptimizer<F>::Optimize(const std::vector<double>& starting_point) {
fncas::reset_internals_singleton();
const size_t dim = starting_point.size();
const fncas::x gradient_helper(dim);
Expand All @@ -194,12 +196,12 @@ OptimizationResult ConjugateGradientOptimizer<F>::Optimize(const std::vector<dou

const auto initial_f_gradf = gi(current_point);
std::vector<double> current_grad = initial_f_gradf.gradient;
std::vector<double> s(current_grad);
fncas::FlipSign(s); // Initial direction = negative gradient.
std::vector<double> s(current_grad); // Direction to search for a minimum.
fncas::FlipSign(s); // Trying first step against the gradient to minimize the function.

for (size_t iteration = 0; iteration < max_steps_; ++iteration) {
// Backtracking line search.
auto new_point = fncas::BackTracking(fi, gi, current_point, s, bt_alpha_, bt_beta_, bt_max_steps_);
const auto new_point = fncas::BackTracking(fi, gi, current_point, s, bt_alpha_, bt_beta_, bt_max_steps_);
const auto new_f_gradf = gi(new_point);

// Calculating direction for the next step.
Expand All @@ -211,7 +213,7 @@ OptimizationResult ConjugateGradientOptimizer<F>::Optimize(const std::vector<dou
current_grad = new_f_gradf.gradient;
current_point = new_point;

// Simple early stopping by the norm of gradient.
// Simple early stopping by the norm of the gradient.
if (std::sqrt(L2Norm(s)) < grad_eps_ && iteration >= min_steps_) {
break;
}
Expand Down
16 changes: 8 additions & 8 deletions test.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -241,10 +241,10 @@ TEST(FNCAS, NaiveGDvsBacktrackingGDOnRosenbrockFunction1000Steps) {
const double x0_err_bt = std::abs(result_bt.point[0] - 1.0);
const double x1_err_n = std::abs(result_naive.point[1] - 1.0);
const double x1_err_bt = std::abs(result_bt.point[1] - 1.0);
EXPECT_TRUE(fncas::IsNormal(x0_err_n));
EXPECT_TRUE(fncas::IsNormal(x1_err_n));
EXPECT_TRUE(fncas::IsNormal(x0_err_bt));
EXPECT_TRUE(fncas::IsNormal(x1_err_bt));
ASSERT_TRUE(fncas::IsNormal(x0_err_n));
ASSERT_TRUE(fncas::IsNormal(x1_err_n));
ASSERT_TRUE(fncas::IsNormal(x0_err_bt));
ASSERT_TRUE(fncas::IsNormal(x1_err_bt));
EXPECT_TRUE(x0_err_bt < x0_err_n);
EXPECT_TRUE(x1_err_bt < x1_err_n);
EXPECT_NEAR(1.0, result_bt.point[0], 1e-6);
Expand All @@ -264,10 +264,10 @@ TEST(FNCAS, ConjugateGDvsBacktrackingGDOnRosenbrockFunction100Steps) {
const double x0_err_bt = std::abs(result_bt.point[0] - 1.0);
const double x1_err_cg = std::abs(result_cg.point[1] - 1.0);
const double x1_err_bt = std::abs(result_bt.point[1] - 1.0);
EXPECT_TRUE(fncas::IsNormal(x0_err_cg));
EXPECT_TRUE(fncas::IsNormal(x1_err_cg));
EXPECT_TRUE(fncas::IsNormal(x0_err_bt));
EXPECT_TRUE(fncas::IsNormal(x1_err_bt));
ASSERT_TRUE(fncas::IsNormal(x0_err_cg));
ASSERT_TRUE(fncas::IsNormal(x1_err_cg));
ASSERT_TRUE(fncas::IsNormal(x0_err_bt));
ASSERT_TRUE(fncas::IsNormal(x1_err_bt));
EXPECT_TRUE(x0_err_cg < x0_err_bt);
EXPECT_TRUE(x1_err_cg < x1_err_bt);
EXPECT_NEAR(1.0, result_cg.point[0], 1e-6);
Expand Down

0 comments on commit 3dac18f

Please sign in to comment.