In [10]:
from math import exp
from scipy import optimize
from scipy.optimize import minimize
from concurrent.futures import ThreadPoolExecutor
import matplotlib.pyplot as plt
import time
import numpy
import itertools
In [11]:
# f(x) = x / (1+x^2), x* = -1
def objfun1(x, xprev, lamb):
return x / (1+x**2) + lamb * abs(x - xprev)**2
# f(x) = -x / e^x, x* = 1
def objfun2(x, xprev, lamb):
return -x / exp(x) + lamb * abs(x - xprev)**2
# f(x) = (x-2)^4 + (x-2)^3 + (x-2)^2 + 1, x* = 2
def objfun3(x, xprev, lamb):
return (x-2)**4 + (x-2)**3 + (x-2)**2 + 1 + lamb * abs(x - xprev)**2
# ff(x) = (x-1)^2 * (x-3)^2 * (x-5)^2 + 1, x1 = 1, x2 = 3 e x3 = 5
def objfun4(x, xprev, lamb):
return (x-1)**2 * (x-3)**2 * (x-5)**2 + 1 + lamb * abs(xprev-x)**2
# f(x) = ((x-1)*(x+1)*(x-2)*(x+2))^2 - 2, x1 = -1, x2 = 1, x3 = -2, x4 = 2
def objfun5(x, xprev, lamb):
return ((x-1)*(x+1)*(x-2)*(x+2))**2 - 2 + lamb * abs(x - xprev)**2
In [70]:
# x, x0, iterações, erro, tempo
XK, X0, IT, ERR, TIME = [], [], [], [], []
def minimalg_single(str_func, it, lambda_str, tolerance, optimum, low, high):
#optimum = float(input("Optimum: "))
#optimum = 2
#tests = int(input("Number of tests: "))
tests = 100
#max_iter = int(input("Max iterations: "))
max_iter = 300
#low = float(input("x0 interval (lower): "))
#high = float(input("x0 interval (upper): "))
objfun = globals()[str_func]
sum_iter = 0
sum_time = 0
sum_x = 0
sum_y = 0
for i in range(1, tests + 1):
x0 = numpy.random.uniform(low, high)
xprev = x0
start_time = time.time()
for k in range(1, max_iter + 1):
lamb = eval(lambda_str[it])
res = minimize(objfun, x0, args=(xprev, lamb), method='BFGS', tol=tolerance[it], options={'disp': False, 'maxiter': 100})
error = abs(xprev-res.x[0])
xprev = res.x[0]
if (error <= tolerance[it]):
break
end_time = time.time() - start_time
XK.append(xprev)
X0.append(x0)
IT.append(k)
ERR.append(error)
TIME.append(end_time)
sum_iter = sum_iter + k
sum_time = sum_time + end_time
sum_x = sum_x + res.x[0]
sum_y = sum_y + res.fun
print('%d & %.5f & $%s$ & %.2f & %.10f & %.10f & %.10f \\\\ [1ex]' % (it+1, tolerance[it], lambda_str[it], sum_iter/tests, sum_time/tests, sum_x/tests, abs((sum_x/tests)-optimum)));
In [38]:
lambda_str = ['3+(-1)**k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k']
tolerance = [1e-1, 1e-2, 1e-3, 1e-4, 1e-5, 1e-6, 1e-7, 1e-8, 1e-9, 1e-10]
# f(x) = x / (1+x^2)
for it in range(0, 10):
minimalg_single("objfun1", it, lambda_str, tolerance, -1, -1, 0.5)
In [68]:
plt.xlabel('Tolerância')
plt.ylabel('Erro')
er = [0.46655473157864590927, 0.10577582452102307631, 0.01400086485775875467, 0.00146590855738426562, 0.00014559136557190655, 0.00001470499862710195, 0.00000136785759341507, 0.00000004398226505220, 0.00000004383897622695, 0.00000005511787470880]
ite = [2.05, 12.21, 33.40, 57.64, 85.30, 110.36, 136.18, 164.72, 180.23, 217.51]
t = [1e-1, 1e-2, 1e-3, 1e-4, 1e-5, 1e-6, 1e-7, 1e-8, 1e-9, 1e-10]
plt.gca().set_xscale('log')
x, y = zip(*sorted(zip(er, t)))
plt.plot(y, x, lw=1.5)
plt.show()
In [40]:
lambda_str = ['1+1/k', '1+1/k', '1+1/k', '1/k', '1/k', '1/k', '1', '1', '1', '2-1/k', '2-1/k', '2-1/k', '3+(-1)**k', '3+(-1)**k', '3+(-1)**k']
tolerance = [1e-3, 1e-4, 1e-5, 1e-3, 1e-4, 1e-5, 1e-3, 1e-4, 1e-5, 1e-3, 1e-4, 1e-5, 1e-3, 1e-4, 1e-5]
# f(x) = x / (1+x^2)
for it in range(0, 15):
minimalg_single("objfun1", it, lambda_str, tolerance, -1, -1, 0.5)
In [71]:
# f(x) = -x / e^x
for it in range(0, 15):
minimalg_single("objfun2", it, lambda_str, tolerance, 1, 1, 3)
In [72]:
# f(x) = (x-2)^4 + (x-2)^3 + (x-2)^2 + 1, x* = 2
for it in range(0, 15):
minimalg_single("objfun3", it, lambda_str, tolerance, 2, -10, 12)
In [122]:
def minimalg_multiple(str_func, lambda_str, it):
optimum = 2
x0 = 5
max_iter = 200
objfun = globals()[str_func]
tolerance = 1e-5
xprev = x0
#lambda_str = input("Lambda function: ")
start_time = time.time()
for k in range(1, max_iter + 1):
lamb = eval(lambda_str[it])
res = minimize(objfun, x0, args=(xprev, lamb), method='BFGS', tol=1e-5, options={'disp': False, 'maxiter': 100})
error = abs(xprev-res.x[0])
xprev = res.x[0]
if (error <= tolerance):
break;
end_time = time.time() - start_time
print('%d & %d & %.1f & $%s$ & %d & %.10f & %.10f & %.10f \\\\ [1ex]' % (it+1, optimum, x0, lambda_str[it], k, end_time, xprev, abs(xprev-optimum)))
In [110]:
lamb = ['1+1/k', '1/k', '1', '2-1/k', '3+(-1)**k']
for i in range(0, 5):
minimalg_multiple("objfun4", lamb, i)
In [123]:
lamb = ['1+1/k', '1/k', '1', '2-1/k', '3+(-1)**k']
for i in range(0, 5):
minimalg_multiple("objfun5", lamb, i)
In [ ]: