黄金分割法

# -*- coding: UTF-8 -*-
from math import *
import numpy as np

############################################
# 黄金分割法(関数の最小値)
#      a,b : 初期区間 a < b
#      eps : 許容誤差
#      val : 間数値
#      ind : 計算状況
#              >= 0 : 正常終了(収束回数)
#              = -1 : 収束せず
#      max : 最大試行回数
#      fun : 関数値を計算する関数の名前
#      return : 結果
#      coded by Y.Suganuma
############################################

def gold(a, b, eps, val, ind, max, fun) :

	x      = 0.0
	tau    = (sqrt(5.0) - 1.0) / 2.0
	count  = 0
	ind[0] = -1
	x1     = b - tau * (b - a)
	x2     = a + tau * (b - a)
	f1     = fun(x1)
	f2     = fun(x2)

	while count < max and ind[0] < 0 :
		count += 1
		if f2 > f1 :
			if (abs(b-a) < eps) and (abs(b-a) < eps*abs(b)) :
				ind[0] = 0
				x      = x1
				val[0] = f1
			else :
				b  = x2
				x2 = x1
				x1 = a + (1.0 - tau) * (b - a)
				f2 = f1
				f1 = fun(x1)
		else :
			if (abs(b-a) < eps) and (abs(b-a) < eps*abs(b)) :
				ind[0] = 0
				x      = x2
				val[0] = f2
				f1     = f2
			else :
				a  = x1
				x1 = x2
				x2 = b - (1.0 - tau) * (b - a)
				f1 = f2
				f2 = fun(x2)

	if ind[0] == 0 :
		ind[0] = count
		fa     = fun(a)
		fb     = fun(b)
		if fb < fa :
			a  = b
			fa = fb
		if fa < f1 :
			x      = a
			val[0] = fa

	return x

----------------------------------

# -*- coding: UTF-8 -*-
import numpy as np
import sys
from math import *
from function import gold

############################################
# 黄金分割法によるy=x^4+3x^3+2x^2+1の最小値
#      coded by Y.Suganuma
############################################

			# 関数値の計算
def snx(x) :
	return x * x * x * x + 3.0 * x * x * x + 2.0 * x * x + 1.0
			# 設定と実行
a   = -2.0
b   = -1.0
eps = 1.0e-10
max = 100
val = np.empty(1, np.float)
ind = np.empty(1, np.int)

x = gold(a, b, eps, val, ind, max, snx)

print("x " + str(x) + " val " + str(val[0]) + " ind " + str(ind[0]))