Broyden's method for solving systems of nonlinear equations.
// classic Broyden, without jacobian bool broyden(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, arma::vec* grad_out, void* opt_data)> opt_objfn, void* opt_data); bool broyden(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, arma::vec* grad_out, void* opt_data)> opt_objfn, void* opt_data, algo_settings_t& settings); // classic Broyden, with jacobian bool broyden(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, void* opt_data)> opt_objfn, void* opt_data, std::function<arma::mat (const arma::vec& vals_inp, void* jacob_data)> jacob_objfn, void* jacob_data); bool broyden(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, void* opt_data)> opt_objfn, void* opt_data, std::function<arma::mat (const arma::vec& vals_inp, void* jacob_data)> jacob_objfn, void* jacob_data, algo_settings_t& settings); // derivative-free method of Li and Fukushima (2000) bool broyden_df(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, void* opt_data)> opt_objfn, void* opt_data); bool broyden_df(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, void* opt_data)> opt_objfn, void* opt_data, algo_settings_t& settings); // derivative-free method with jacobian bool broyden_df(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, void* opt_data)> opt_objfn, void* opt_data, std::function<arma::mat (const arma::vec& vals_inp, void* jacob_data)> jacob_objfn, void* jacob_data); bool broyden_df(arma::vec& init_out_vals, std::function<arma::vec (const arma::vec& vals_inp, void* opt_data)> opt_objfn, void* opt_data, std::function<arma::mat (const arma::vec& vals_inp, void* jacob_data)> jacob_objfn, void* jacob_data, algo_settings_t& settings);
Function arguments:
init_out_vals
a column vector of initial values; will contain the final values.opt_objfn
the function to be zeroed-out, taking two arguments:
vals_inp
a vector of inputs; andopt_data
additional parameters passed to the function.opt_data
additional parameters passed to the function.jacob_objfn
a function producing the jacobian matrix, taking two arguments:
vals_inp
a vector of inputs; andjacob_data
additional parameters passed to the function.jacob_data
additional parameters passed to the function.settings
parameters controlling the optimization routine; see below.Optimization control parameters:
double err_tol
the value controlling how small $\| f \|$ should be before 'convergence' is declared.int iter_max
the maximum number of iterations/updates before the algorithm exits.There are two forms of Broyden's method available in OptimLib.
Let $x^{(i)}$ denote the function input values at stage $i$ of the algorithm. The Broyden updating rule is as follows.
Let
$$s^{(i)} := x^{(i)} - x^{(i-1)}$$ $$y^{(i)} := F(x^{(i)}) - F(x^{(i-1)})$$The algorithm stops when $\| F \|$ is less than err_tol
, or the total number of iterations exceeds iter_max
.
Let $x^{(i)}$ denote the function input values at stage $i$ of the algorithm. The updating rule of Li and Fukushima (2000) is as follows.
Let
$$s^{(i)} := x^{(i)} - x^{(i-1)}$$ $$y^{(i)} := F(x^{(i)}) - F(x^{(i-1)})$$The algorithm stops when $\| F \|$ is less than err_tol
, or the total number of iterations exceeds iter_max
.
This example is taken from Matlab's website.
$$F(x) = \begin{bmatrix} \exp(-\exp(-(x_1+x_2))) - x_2(1+x_1^2) \\ x_1\cos(x_2) + x_2\sin(x_1) - 0.5 \end{bmatrix}$$The solution is $(0.3532,0.6061)$.
Code to solve this problem, with and without the jacobian, is given below.
#include "optim.hpp" arma::vec zeros_test_objfn_1(const arma::vec& vals_inp, void* opt_data) { double x_1 = vals_inp(0); double x_2 = vals_inp(1); arma::vec ret(2); ret(0) = std::exp(-std::exp(-(x_1+x_2))) - x_2*(1 + std::pow(x_1,2)); ret(1) = x_1*std::cos(x_2) + x_2*std::sin(x_1) - 0.5; // return ret; } arma::mat zeros_test_jacob_1(const arma::vec& vals_inp, void* opt_data) { double x_1 = vals_inp(0); double x_2 = vals_inp(1); arma::mat ret(2,2); ret(0,0) = std::exp(-std::exp(-(x_1+x_2))-(x_1+x_2)) - 2*x_1*x_1; ret(0,1) = std::exp(-std::exp(-(x_1+x_2))-(x_1+x_2)) - x_1*x_1 - 1.0; ret(1,0) = std::cos(x_2) + x_2*std::cos(x_1); ret(1,1) = -x_1*std::sin(x_2) + std::cos(x_1); // return ret; } int main() { // // Broyden Derivative-free Method of Li and Fukushima (2000) arma::vec x_1 = arma::zeros(2,1); bool success_1 = optim::broyden_df(x_1,zeros_test_objfn_1,nullptr); if (success_1) { std::cout << "broyden_df: test_1 completed successfully." << std::endl; } else { std::cout << "broyden_df: test_1 completed unsuccessfully." << std::endl; } arma::cout << "broyden_df: solution to test_1:\n" << x_1 << arma::endl; // // Derivative-free Method of Li and Fukushima (2000) using the jacobian x_1 = arma::zeros(2,1); success_1 = optim::broyden_df(x_1,zeros_test_objfn_1,nullptr,zeros_test_jacob_1,nullptr); if (success_1) { std::cout << "broyden_df w jacobian: test_1 completed successfully." << std::endl; } else { std::cout << "broyden_df w jacobian: test_1 completed unsuccessfully." << std::endl; } arma::cout << "broyden_df w jacobian: solution to test_1:\n" << x_1 << arma::endl; // // standard Broyden method x_1 = arma::zeros(2,1); success_1 = optim::broyden(x_1,zeros_test_objfn_1,nullptr); if (success_1) { std::cout << "broyden: test_1 completed successfully." << std::endl; } else { std::cout << "broyden: test_1 completed unsuccessfully." << std::endl; } arma::cout << "broyden: solution to test_1:\n" << x_1 << arma::endl; // // standard Broyden method using the jacobian x_1 = arma::zeros(2,1); success_1 = optim::broyden(x_1,zeros_test_objfn_1,nullptr,zeros_test_jacob_1,nullptr); if (success_1) { std::cout << "broyden w jacobian: test_1 completed successfully." << std::endl; } else { std::cout << "broyden w jacobian: test_1 completed unsuccessfully." << std::endl; } arma::cout << "broyden w jacobian: solution to test_1:\n" << x_1 << arma::endl; return 0; }