Broyden's method for solving systems of nonlinear equations.


Definition and Syntax

// 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; and
    • opt_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; and
    • jacob_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.

Details

There are two forms of Broyden's method available in OptimLib.


Basic Broyden's method

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)})$$
  • Update the inverse Jacobian matrix.
  • $$B^{(i)} = B^{(i-1)} + \frac{(s^{(i)} - B^{(i-1)}y^{(i)})y^{(i)}{'}}{y^{(i)} \cdot y^{(i)}}$$
  • Compute the descent direction.
  • $$d^{(i)} = - B^{(i)} F(x^{(i)})$$
  • Update the input values.
  • $$x^{(i+1)} = x^{(i)} + d^{(i)}$$

The algorithm stops when $\| F \|$ is less than err_tol, or the total number of iterations exceeds iter_max.


Derivative-free Method of Li and Fukushima (2000)

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)})$$
  • Update the inverse Jacobian matrix.
  • $$B^{(i)} = B^{(i-1)} + \frac{(s^{(i)} - B^{(i-1)}y^{(i)})y^{(i)}{'}}{y^{(i)} \cdot y^{(i)}}$$
  • Compute the descent direction.
  • $$d^{(i)} = - B^{(i)} F(x^{(i)})$$
  • Quasi line search step.
  • If $$\| F(x^{(i)} + d^{(i)}) \| \leq \rho \| F(x^{(i)}) \| - \sigma_2 \| d^{(i)} \|^2$$ then:
    • set $\lambda_i = 1.0$
    otherwise:
    • let $\eta_i = 1 / i^2$, and find the smallest $k$ for which $$\| F(x^{(i)} + \lambda_i d^{(i)}) \| \leq \| F(x^{(i)}) \| - \sigma_1 \| d^{(i)} \|^2 + \eta_i \| F(x^{(i)}) \|$$ holds with $\lambda_i = \beta^k$.

      Control parameter values: $\rho = 0.9$, $\sigma_1 = 0.001$, and $\sigma_2 = 0.001$.
  • Update the input values.
  • $$x^{(i+1)} = x^{(i)} + \lambda_i d^{(i)}$$

The algorithm stops when $\| F \|$ is less than err_tol, or the total number of iterations exceeds iter_max.


Examples

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;
}