実数探索三種類解説
自分がよく使う実数上の探索アルゴリズム「三分探索」「黄金分割探索」「二分探索」のメモです。
三分探索
三分探索は凸関数の極値を求めるために使うアルゴリズムです。このアルゴリズムは関数が微分不可能な場合にも使えます。やり方は探索領域を三分割し、真ん中二本の境界のうちどちらの値が大きいかを調べ、それとは逆の方の境界を新たな探索領域の端にするというものです。一回のイテレーションで関数の計算を2回行い、探索領域が2/3になります。
double search(double left, double right) { for (int loop = 0; loop < maxLoop; ++loop){ if (f((left * 2 + right) / 3) > f((left + right * 2) / 3)){ right = (left + right * 2) / 3; } else { left = (left * 2 + right) / 3; } } return (right + left) * 0.5; }
黄金分割探索
黄金分割探索は三分探索と同様に凸関数の極値を求めるために使うアルゴリズムです。領域の分割を黄金比によって再帰的に行うため、前回のイテレーションでの値を再利用することができ、高速化することができます。1回のイテレーションで関数の計算を1回行い、探索領域が約0.618倍になります。
static const double RATIO = (1 + sqrt(5.0)) * 0.5; double search(double left, double rithg){ double leftValue = f((left * RATIO + right) / (1 + RATIO)); double rightValue = f((left + right * RATIO) / (1 + RATIO)); for (int loop = 0; loop < maxLoop; ++loop){ if (leftValue > rightValue){ right = (left + right * RATIO) / (1 + RATIO); rightValue = leftValue; leftValue = f((left * RATIO + right) / (1 + RATIO)); } else { left = (left * RATIO + right) / (1 + RATIO); leftValue = rightValue; rightValue = f((left + right * RATIO) / (1 + RATIO)); } } if (leftValue > rightValue){ right = (left + right * RATIO) / (1 + RATIO); } else { left = (left * RATIO + right) / (1 + RATIO); } return (right + left) * 0.5; }
二分探索
二分探索(実数上)は単調増加(減少)関数がx=kと交わる点を求めるために使うアルゴリズムです。下記の例ではk=0としています。微分可能な凸関数の極大を求める場合は、凸関数を微分した関数をf()としてください。1回のイテレーションで関数の計算を1回行い、探索領域が約0.5倍になります。
double search(double left, double right) { for (int loop = 0; loop < maxLoop; ++loop){ if (f((left + right) / 2) > 0){ right = (left + right) / 2; } else { left = (left + right) / 2; } } return (left + right) / 2; }
注意点
よくループの終了条件を
while (left - right > EPS)
としているコードを見かけますが、double(64bit浮動小数)の有効数字は53bit分、十進数で15桁分の為、値が大きいときに終了しなくなってしまう場合があります。これを避けるため、必ずループ回数を固定して探索を行うことを強く推奨いたします。64bit浮動小数を過不足無く収束させるための基準となる回数は以下の通りです。
探索アルゴリズム | ループ回数の目安 | 導出方法 |
三分探索 | 86回 | (2/3)^{n}<10^{-15} |
黄金分割探索 | 72回 | ((-1+√5)/2)^{n}<10^{-15} |
二分探索 | 50回 | (1/2)^{n}<10^{-15} |