<?php
/*********************************************/
/* 黄金分割法によるy=x^4+3x^3+2x^2+1の最小値 */
/* coded by Y.Suganuma */
/*********************************************/
$a = -2.0;
$b = -1.0;
$eps = 1.0e-10;
$max = 100;
$x = gold($a, $b, $eps, $val, $ind, $max, "snx");
printf("x %f val %f ind %d\n", $x, $val, $ind);
/****************/
/* 関数値の計算 */
/****************/
function snx($x)
{
return $x * $x * $x * $x + 3.0 * $x * $x * $x + 2.0 * $x * $x + 1.0;
}
/******************************************/
/* 黄金分割法(関数の最小値) */
/* a,b : 初期区間 a < b */
/* eps : 許容誤差 */
/* val : 間数値 */
/* ind : 計算状況 */
/* >= 0 : 正常終了(収束回数) */
/* = -1 : 収束せず */
/* max : 最大試行回数 */
/* fun : 関数値を計算する関数の名前 */
/* return : 結果 */
/******************************************/
function gold($a, $b, $eps, &$val, &$ind, $max, $fun)
{
$x = 0.0;
$tau = (sqrt(5.0) - 1.0) / 2.0;
$count = 0;
$ind = -1;
$x1 = $b - $tau * ($b - $a);
$x2 = $a + $tau * ($b - $a);
$f1 = $fun($x1);
$f2 = $fun($x2);
while ($count < $max && $ind < 0) {
$count += 1;
if ($f2 > $f1) {
if (abs($b-$a) < $eps && abs($b-$a) < $eps*abs($b)) {
$ind = 0;
$x = $x1;
$val = $f1;
}
else {
$b = $x2;
$x2 = $x1;
$x1 = $a + (1.0 - $tau) * ($b - $a);
$f2 = $f1;
$f1 = $fun($x1);
}
}
else {
if (abs($b-$a) < $eps && abs($b-$a) < $eps*abs($b)) {
$ind = 0;
$x = $x2;
$val = $f2;
$f1 = $f2;
}
else {
$a = $x1;
$x1 = $x2;
$x2 = $b - (1.0 - $tau) * ($b - $a);
$f1 = $f2;
$f2 = $fun($x2);
}
}
}
if ($ind == 0) {
$ind = $count;
$fa = $fun($a);
$fb = $fun($b);
if ($fb < $fa) {
$a = $b;
$fa = $fb;
}
if ($fa < $f1) {
$x = $a;
$val = $fa;
}
}
return $x;
}
?>