# -*- coding: utf-8 -*-
# mathtoolspy
# -----------
# A fast, efficient Python library for mathematically operations, like
# integration, solver, distributions and other useful functions.
#
# Author: sonntagsgesicht, based on a fork of Deutsche Postbank [pbrisk]
# Version: 0.3, copyright Wednesday, 18 September 2019
# Website: https://github.com/sonntagsgesicht/mathtoolspy
# License: Apache License 2.0 (see LICENSE file)
from mathtoolspy.utils.mathconst import DOUBLE_TOL, GOLD
from mathtoolspy.utils.math_fcts import float_equal, abs_sign
[docs]def shift(a, b, c):
''' The shift function returns the values as tuple.'''
return a, b, c
# devFct, tInitialValue, tolerance, Optimizer.MAX_NUMBER_OF_FUNCTION_CALLS
[docs]def minimize_algorithm_1dim_brent(fct, _a, _b, _c, tolerance=DOUBLE_TOL):
'''
Finds the minimum of the given function f. The arguments are the given function f, and given a bracketing triplet of abscissas A, B, C
(such that B is between A and C, and f(B) is less than both f(A) and f(C)) and the Tolerance.
This routine isolates the minimum to a fractional precision of about tol using
Brent's method. The abscissa of of the minimum is returned as xmin, and the minimum
value is returned as brent, the returned function value.
'''
''' ZEPS is a small number that protects against trying to achieve fractional accuracy
for a minimum that happens to be exactly zero. '''
ZEPS = 1.0e-10
a = _a if _a < _c else _c
b = _a if _a > _c else _c
if not a < _b < b:
raise RuntimeError("Value %0.4f does not embrace %0.4f and %0.4f for bracketing." % (_b, a, b))
x = w = v = _b
fv = fw = fx = fct(x)
tol1 = tolerance
d = e = 0.0
e_temp = fu = u = xm = 0.0
iterations = 0
while (True):
xm = 0.5 * (a + b)
tol1 = tolerance * abs(x) + ZEPS
tol2 = 2.0 * tol1;
if abs(x - xm) <= tol2 - 0.5 * (b - a):
return (x, fx)
if abs(e) > tol1:
r = (x - w) * (fx - fv)
q = (x - v) * (fx - fw)
p = (x - v) * q - (x - w) * r
q = 2.0 * (q - r)
if q > 0.0:
p = -p
q = abs(q)
e_temp = e
e = d
if abs(p) >= abs(0.5 * q * e_temp) or p <= q * (a - x) or p >= q * (b - x):
e = a - x if x >= xm else b - x
d = GOLD * e
else:
d = p / q
u = x + d
if u - a < tol2 or b - u < tol2:
d = abs_sign(tol1, xm - x)
else:
e = a - x if x >= xm else b - x
d = GOLD * e
u = x + d if abs(d) >= tol1 else x + abs_sign(tol1, d)
fu = fct(u);
if fu <= fx:
if u >= x:
a = x
else:
b = x
v, w, x = shift(w, x, u)
fv, fw, fx = shift(fw, fx, fu)
else:
if u < x:
a = u
else:
b = u
if fu <= fw or float_equal(w, x):
v = w
w = u
fv = fw
fw = fu
elif float_equal(fu, fv) or float_equal(v, x) or float_equal(v, w):
v = u
fv = fu
iterations = iterations + 1
if iterations > 10000:
return (None, None)