2次関数の区間内での極大値

任意の2次関数f(x)について、[x1, x2]における極大値を探すプログラムをC言語で実装した。

#include <stdio.h>

double some_quadratic_function(double x) {
    return -2.2 * x * x + 4.686 * x + 5.467;
}

double quadratic_max(double (*f)(double), double x1, double x2) {
    double epsilon = 1e-7;
    double xa, xb, xc;

    // 探索
    xa = x1 < x2 ? x1 : x2;
    xb = x1 < x2 ? x2 : x1;
    while (xb - xa > epsilon) {
        xc = (xa + xb) / 2.0;
        if ( f(xa) < f(xb) )
            xa = xc;
        else
            xb = xc;
    }

    return f(x1) < f(x2) ? ( f(x2) < f(xc) ? f(xc) : f(x2) ) : f(x1);
}

int main(void) {
    printf("%lf\n", quadratic_max(some_quadratic_function, -10.0, 10.0));
    return 0;
}

まず、任意の2次関数を適当に宣言している。2次関数でなくてもよいが、それっぽい関数というだけではうまく探せるとは限らない(後述)。

quadratic_max()がこのプログラムのキモとなる関数だ。まずは第一引数を見よう。これは、double型の引数をひとつとり、double型の値を返す関数のポインタにfという名前をつける、という意味である。C言語の関数ポインタはこんな書き方で表現する。
次に注目するのは、極大値を求める探索だ。といっても2次関数について調べる設計であるため、たいした探索アルゴリズムではない。xa, xbは探索範囲の上下の値である。xbxaより小さくないことがつねに保証されている。そして、f(xa)f(xb)を比較し、f(xb)のほうが大きければ次回の探索をxcxbで、そうでなければxaxcで行う。ここでxcxa, xbの中間点である。探索を繰り返すと範囲が狭まるので、探索範囲がepsilonより小さくなったら探索を終える。最後に、今求めた点xcでの値f(xc)や元の探索範囲の境界の値f(x1), f(x2)を比べ、最も大きかった値を返す。

ユーザーは、調べたい関数の仮引数のラベル(f(x)のx, C(t)のt)を考慮せずに、関数の微分形を知らずにこのコードを実装できる。これがこのプログラムの最大のメリットである。

また極小値を知りたいときは、quadratic_maxを書き換えてもよいが、minus_f()を新たに定義し、次のようにすることで、極大値問題に置きかえることができる。

double minus_f(double x) { return - some_quadratic_function(x); }

int main(void) {
    printf("%lf\n", quadratic_max(minus_f, -10.0, 10.0));
    return 0;
}

前に述べたように、このアルゴリズムを2次関数に対して用いる限りは問題ない。しかし、2次関数と似た「探索範囲内で連続であり、唯一の極大値を持つ関数」であっても、適用できない関数が存在する。たとえばx\exp(-x)について[-1, 9]で探索すると、間違いなく失敗する(下図)。それはアルゴリズムの中で、2次関数の持つ対称性を暗黙のうちに利用しているためである。