読者です 読者をやめる 読者になる 読者になる

競技プログラミングをするんだよ

ICPC国内予選突破を目標に一日一問題以上解いていきます。

凸関数(単峰関数)の極値の探索について

実装メモ

単峰関数の極値の探索についてです.

あんまり数学的な記述はしたくない(プロに殺されかねない)ので,単峰関数の説明は省きます.
凸関数のような極値を一つだけもつお山の形をした関数です.

競技プログラミングにおいて,こういった問題を解くときに使われるアルゴリズム

  • 三分探索
  • 黄金分割探索 (実は使ったことがない)

です.
三分探索は解を含む区間を三等分して区間を狭めていく方法です.
以下のようなコードで実装できます.(f(x)の最小値をとるxを求める)

double lb=-INF,ub=INF;
rep(i, 100) {
	double ml = (lb + lb + ub) / 3;
	double mu = (lb + ub + ub) / 3;
	if (f(ml) < f(mu))ub = mu;
	else lb = ml;
}

で,なんですが
三分探索は一回の操作で区間の長さが2/3になります.
計算量はlogになるのですが,二分探索などであらわれるlogよりも底の値が小さくなるので精度を求めるために必要なループの回数が大体2倍ぐらいになります.

これでも十分実用的なのですが,極値の探索を底が純粋な2のlogで探索する方法があるので,紹介したいと思います.
一つの方法は微分して二分探索なのですが,数学の苦手な私には使いこなせないのでお勧めしません.

紹介する方法は,山登り法の歩幅を半分にしていくようなイメージです.
大学のラボの同期と問題を解いていたときに偶然編み出した方法なのですが,114514回既出かもしれないです(未検証)

アルゴリズム

以下のようなアルゴリズムで実現します.

  1. x(初期値は何でもよい) ,w=inf
  2. x=argmin_{x \in \left\{x-w,x,x+w \right\}}f(x) , w=w/2
  3. 2の繰り返し

証明

上手く解が発見できる理由について証明します.間違ってたら優しく指摘してください

f(x)の極小値をとるxx_a,現在のxx_tとします.


常に区間\left[x_t-2w,x_t+2w \right]の中に解x_aが含まれていることを示します.
初期状態において,w=infとしているので,この条件を必ず満たします.

  • x_t \le x_a \le x_t+wのとき

遷移後のx'_tについて,凸性よりf(x_t-w) > f(x_t)が成りたつので,x'_t=x_t-wになることはありません.
このとき,x_t,x_t+wのどちらに遷移させても,|x'_t - x_a| \le wが成り立ちます.
したがって,遷移後についても|x'_t-x_a|\le 2w'(w'=\frac{w}{2}) が成り立ちます.

  • x_t \ge x_a \ge x_t-wのとき

↑のを反対にするだけです.

  • x_t+w \le x_a \le x_t+2wのとき

このとき,必ずx'_t=x_t+wとなります.
そのため,この場合の遷移後も|x'_t-x_a|\le 2w'(w'=\frac{w}{2}) が成り立ちます.

  • x_t-w \ge x_a \ge x_t-2wのとき

↑のを反対にするだけです.

よって,解が幅4w区間内に含まれ,遷移を繰り返すごとに,区間のサイズが半分になっていくため,区間が2/3になる三分探索よりも早く収束します.
また,三分探索では毎ラウンド二つの点について値を計算するのですが,この解法ではifelseを用いることによって1回のラウンドで一つの値の計算のみで済んだりします(最悪は二つ)

実装例

double x=s,w=INF;
double y=f(s),nx,ny;
rep(i, 100) {
		if ((ny=f(nx=x+w)) < y) y = ny, x = nx; 
	else	if ((ny=f(nx=x-w)) < y) y = ny, x = nx;
	w /= 2;
}