Newton's method in C ++, Python, Go (for understanding function objects)

Overview

Recently, I often see implementations that pass function objects as arguments of callback functions in C ++, and I'm still confused, so I want to implement something myself and fix it in my head. This time, I chose Newton's method as a simple calculation formula. Implement Newton's method in C ++, Python, Go using function objects.

Target audience

--C ++ programmer --Python programmer --Go programmer --A person who understands the basics of programming

Function object

A function object is any object that has a function call operator defined. From cppreference

Articles like the one below also explain function objects and callback functions. Inline expansion of function pointers and function objects Series introducing small C ++ technology [Small C ++ 9 times] # 3

Details

Newton's method

Newton's method is a method that can find a solution even with a non-linear equation by numerical calculation. An algorithm that finds the tangent line from the initial value X1 and performs iterative calculation while updating X at the intersection of the tangent line and the X axis.

Image https://linky-juku.com/linear-programming/

If it is a linear equation, the solution of the equation can be answered immediately by using the proportionality constant. With the formula below, the solution is immediately found to be $ x =-\ frac {1} {2} $.

f(x) = 2x + 1

So what about the formula below?

f(x) = x^3 - 5x + 1

If the equation can be transformed well and can be made like $ (x-a) (x-b) (x-c) $, a solution can be obtained even with a nonlinear equation. This is not possible with actual nonlinear equations (eg, geothermal diffusion equations and turbulence calculations). At such times, there is Newton's method to force iterative calculations to find solutions. The equation of Newton's method is as follows.

x_{n+1}=x_n−\frac{f(x_n)}{f′(x_n)}

I think that Newton's method is often written in programming classes at university. It's not a difficult formula, so it's easy to get started.

In this article as well, I chose the simple Newton's method to learn while writing function objects myself.

Let's implementation

Basic interface

1. 1. Create function objects for f (x) and df (x).
2. Pass the function object created in the Newton function, the initial value of X, the maximum number of trials, and the convergence condition as arguments.

• This time, when the initial value is x = 0, 2, 3, the maximum number of trials is 1000,

3. 3. Implement the following polynomial with the Newton function. It is implemented in an infinite loop, and the calculation ends when the difference between $ x_ {n + 1} $ and $ x_n $ becomes smaller than $ 1.0 × 10 ^ {-10} $. Also, the calculation ends even when the maximum number of trials is reached.

x_{n+1}=x_n−\frac{f(x_n)}{f′(x_n)}

4. Receives the return value of the Newton function and displays the result.

C++

1. 1. Create function objects for f (x) and df (x).

In C ++ we used std :: pow to implement the cube of X. In this case, the formula is not as it is, and I think it is not readable.

Main.cpp

double fx(double x)
{
    return std::pow(x, 3.0) - 5 * x + 1;
}

double df(double x)
{
    return 3 * std::pow(x, 2) - 5;
}

2. Pass the created function object, initial value of X, maximum number of trials, and convergence condition as arguments to the Newton function.

Main.cpp

    Calc calc;
    std::function<double(std::function<double(double)>,
                         std::function<double(double)>,
                         double,
                         int,
                         double)> newton = calc;

    std::cout << "Newton Object : " <<newton(fx, df, 2, Calc::MAX_ITER, Calc::EPC) << std::endl;
    std::cout << "Newton Object : " <<newton(fx, df, 0, Calc::MAX_ITER, Calc::EPC) << std::endl;
    std::cout << "Newton Object : " <<newton(fx, df, 3, Calc::MAX_ITER, Calc::EPC) << std::endl;

Main.cpp

    std::function<double(std::function<double(double)>,
                         std::function<double(double)>,
                         double,
                         int,
                         double)> newton_ptr = Newton_Main;

    std::cout << "Newton Ptr : " << newton_ptr(fx, df, 2, Calc::MAX_ITER, Calc::EPC) << std::endl;
    std::cout << "Newton Ptr : " << newton_ptr(fx, df, 0, Calc::MAX_ITER, Calc::EPC) << std::endl;
    std::cout << "Newton Ptr : " << newton_ptr(fx, df, 3, Calc::MAX_ITER, Calc::EPC) << std::endl;

3. 3. Implement the following polynomial with the Newton function. Implement in an infinite loop so that the calculation ends with the convergence condition and the maximum number of trials.

By defining the () operator in the Calc class,

Calc.cpp

double Calc::operator()(std::function<double(double)>f, std::function<double(double)>df, double x0, int max_iter, double epc)
{
    double x = x0;
    int iter = 0;
    while(1)
    {
        xNew_ = x - f(x)/df(x);
        if (std::abs(x - xNew_) < epc)
        {
            break;
        }
        x = xNew_;
        iter ++;
        if (iter == max_iter)
        {
            break;
        }
    }
    return xNew_;
}

Main.cpp

double Newton_Main(std::function<double(double)>f, std::function<double(double)>df, double x0, int max_iter, double epc)
{
    double xNew = 0;
    double x = x0;
    int iter = 0;
    while(1)
    {
        xNew = x - f(x)/df(x);
        if (std::abs(x - xNew) < epc)
        {
            break;
        }
        x = xNew;
        iter ++;
        if (iter == max_iter)
        {
            break;
        }
    }
    return xNew;
}

4. Receives the return value of the Newton function and displays the result.

Python

1. 1. Create function objects for f (x) and df (x).

main.py

def f(x):
    return x**3 - 5*x + 1

def df(x):
    return 3*x**2 - 5

2. Pass the created function object, initial value of X, maximum number of trials, and convergence condition as arguments to the Newton function.

In 1, only the definition was made in main.py. The defined function is passed as a function object below.

main.py

    newton = formula.newton_func(f, df)

The maximum number of trials and convergence condition were not passed as arguments, but made into member variables of the Calc class.

calc.py

    def __init__(self):
        self.eps = 1e-10
        self.max_iter = 1000

3. 3. Implement the following polynomial with the Newton function. Implement in an infinite loop so that the calculation ends with the convergence condition and the maximum number of trials.

calc.py

    def newton_func(self, f, df):
        def newton(x0):
            x = x0
            iter = 0
            while True:
                x_new = x - f(x)/df(x)
                if abs(x-x_new) < self.eps:
                    break
                x = x_new
                iter += 1
                if iter == self.max_iter:
                    break
            return x_new
        return newton

4. Receives the return value of the Newton function and displays the result.

Go

1. 1. Create function objects for f (x) and df (x).

I couldn't implement cube with x ** 3 like Python.

main.go

// Calc Newton Calculaion struct
type Calc struct{}

// fx f(x) formula
func (calc Calc) fx(x float64) float64 {
	return x*x*x - 5*x + 1
}

// df differentiated f(x)
func (calc Calc) df(x float64) float64 {
	return 3*x*x - 5
}

As shown below, the method of returning a function object in the main function and passing it as an argument to the Newton function could not be implemented well. It is more intuitive and easy to understand if it is defined and implemented as above. Please let me know if you are kind m (_ _) m

main.go

type Calc struct{}

func (calc Calc) fx(x float64) func() float64 {
	return func() float64 {
		return x*x*x - 5*x + 1
	}
}

func (calc Calc) df(x float64) func() float64 {
	return func() float64 {
		return 3*x*x - 5
	}
}

2. Pass the created function object, initial value of X, maximum number of trials, and convergence condition as arguments to the Newton function.

main.go

func main() {
	calc := Calc{}

	var ans float64

	ans = calc.Newton_main(calc.fx, calc.df, 2, 1000, 1e-10)
	fmt.Println("The answer is : ", ans)

	ans = calc.Newton_main(calc.fx, calc.df, 0, 1000, 1e-10)
	fmt.Println("The answer is : ", ans)

	ans = calc.Newton_main(calc.fx, calc.df, 3, 1000, 1e-10)
	fmt.Println("The answer is : ", ans)
}

3. 3. Implement the following polynomial with the Newton function. Implement in an infinite loop so that the calculation ends with the convergence condition and the maximum number of trials.

main.go

// Newton_main Calculation for Newton method
func (calc Calc) Newton_main(fx func(float64) float64, df func(float64) float64, x0 float64, maxIter int, epc float64) float64 {
	var xNew float64
	var x = x0
	var iter int
	for {
		xNew = x - fx(x)/df(x)
		if math.Abs(x-xNew) < epc {
			break
		}
		x = xNew
		iter++
		if iter == maxIter {
			break
		}
	}
	return xNew
}

4. Receives the return value of the Newton function and displays the result.

Summary

I think it is essential for beginners to escape by mastering function objects and closures. Let's learn while using lambdas and closures for answers such as AtCoder. Thank you for reading to the end. The quality of the code is not good, so please comment if you like.

I will add std :: bind later.