SeDuMi の eigK.m がバグってるっぽい？

ここ半年ぐらい、SeDuMi の eigK を実行しているときに良くわからないエラーで止まっていたが、どうもバグっぽい様子。

具体的には、
https://github.com/sqlp/sedumi/blob/master/eigK.m
の７７行目がバグのようで、

	74: li = 0;
	75: xi = nf;
	76: lab = zeros(N,1);
	77: li(li+1:li+nl) = x(xi+1:xi+nl);
	78: xi = xi + nl;

となっていて、74行目でスカラーと定義されている li に77行目で配列代入している。

おそらく、li(li+1:li+nl) ではなくて、lab(li+1:lib+nl) なのではないかと思っていて、自分のところではエラーが消えてなくなったけど、もう少し確認してみたほうがいいのかもしれない。

論文読み：Decision Diagram Decomposition for Quadratically Constrained Binary Optimization

Optimization Online に載っていた以下の論文を読んでみた。

Decision Diagram Decomposition for Quadratically Constrained Binary Optimization
http://www.optimization-online.org/DB_HTML/2018/10/6837.html

扱っている問題は、目的関数が線形で、制約は2次式、あと各変数は0-1変数という問題。
これを MILP に変形して解いている。
この変形に Diagram を使って、うまいことやっている。

基本的には変数の数が増えているので、計算時間が短縮できるかどうかは入力問題しだいだけど、Diagram を使っているところが、なかなかに面白い工夫だと思う。

京都で制御と最適化のワークショップに参加してきた

京都で行われた制御と最適化のワークショップに参加してきた。

SDPについては主要な応用の一つに制御理論があるので、制御分野と数理最適化分野は比較的近い分野にあるわけだけど、一つのワークショップになっているのは、あんまりないかもしれない、と思う。

今回は制御分野の先生からの発表やコメントなどもいろいろと聞くことができて、いつもとは違う視点だったりで参考になるところが多々あった。こういう風に「近い分野だけど、いつもは一緒にやっていない」というところがワークショップを行うのは、いい刺激になるな、と思った。


あと、ワークショップが終わった後に話をしていたところで話題の一つになっていたのが「数理最適化の分野は、一般化しすぎる」ということだった。
分かりやすく錐最適化の言葉で例を作るなら「LPでできたことをSDPとか対称錐とかに一般化する」ということで、つまり「むしろ計算対象を特化することで、より効率的なアルゴリズムを作れるはずだ」ということだった。

数理最適化分野で一般化が良く行われる理由は、個人的には「どういう応用があるかわかっているから」ということだと思っている。つまり、「この問題を解くことに需要はあるのだろうか？」ということを余り真剣に考える必要がない。逆に計算対象を特化すると「この狭い範囲しか計算対象にしないアルゴリズムに、本当に需要はあるんだろうか」ということをキチンと議論しないとならない。

一般化と特化は、どちらかが正解で残りが間違い、という二者択一のものではなくて、両方とも正解であると個人的には思う。どっちに行ってもオモチロイ研究を作り出せることができる。

こういったことを考える視点をもらえた、というだけでも、今回のワークショップは行ってきてよかったなぁ、と思う。
またいつか参加してみたいものだ。

xkeymacs のＯＫボタン

Windows の標準ショートカットキーは使いにくいので（特にカーソル移動をショートカットでできないのが痛い）、xkeymacs を使っている。
ただ、xkeymacs は設定画面が昔のままなので、高い解像度だと「ＯＫ」ボタンがウィンドウからはみ出してしまい、そのままだと押せない。

いろいろと調べてみたら、なんとウィンドウの中に「ＯＫ」ボタンが残っていた。(https://www52.atwiki.jp/wbdev/pages/20.html)
この「ＯＫ」ボタンを発見した人、すごすぎる。いったい、どういう技術を使ったんだろうか、気になる。

ちなみに、「ＯＫ」ボタンを押す変わりに Enter キーを押しても設定できる。

OPTAでもらった質問を考えてみた

ずっと前に OPTA で発表させてもらったときに
「定式化としては線形の目的関数と凸２次制約を含む MI-SOCP だけど、ポートフォリオの計算のように凸２次を目的関数にして線形制約を入れる MI-QP にすれば、もっと短時間の計算でできるのでは？」
というコメントを頂いた。

せっかくいただいたコメントを無駄にするのももったいない、ということで、試行錯誤を繰り返してみること４か月。結論をまず書くと
「短時間にならない。。。」
ということだった。

MI-SOCP で解いた時も、分枝限定法の上界と下界のギャップが小さくできるまでに必要なノードの数が非常に増えてしまったが、MI-QP で行ってもノード数がやはり増えてしまった。特にギャップが2%を下回ると、そこからの改善がどうにも得られなくてノード数の増加を招いている。

今回数値計算している問題は、凸２次のところに現れる係数行列が特殊な構造を持っているので、そのあたりが解くのを難しくしているのかもしれない。
また、いろいろと試行錯誤してみることにする。

Julia は関数に入れると高速化できる

最近になってシミュレーテッド・アニーリングのプログラムを作って C と Julia で比較したところ、C は 1.2 秒ぐらいなのに Julia だと 30 秒ぐらいかかっていて、「Julia、遅いなぁ～」と思っていた。

ところが、
https://qiita.com/ceptree/items/b5cfca180e85e61c42a8
によると、「関数に入れないと Julia では速くならない」ということなので、さっそくプログラムの最初に

function test()

と入れて、最後に
end # end of test()
test()

として
@time include("test.jl")
としてみたら、ほぼ C 言語と同じ 1.2 秒台まで短くなった。

何回も繰り返して実行してみると、C言語のほうが速かったり Julia のほうが速かったり、という感じで、言語の差というよりも他プロセスの動作状況のほうが影響が大きい感じだった。

Julia 1.0 が出てた

Julia 1.0 が８月上旬にリリースされていた。
早速インストールしてみたけど、まだまだ使える状態ではないので、Julia 0.6 か 0.7 のほうが良さそう。

Major version up ということで、言語仕様が変わっていて、周りのパッケージがまだ追い付いていない状況であって、さっきやってみた感じでは CPLEX のパッケージなどは使えない。
(いままで is_linux() とかだった関数が Sys.islinux() などに変更になっていたりする。）

あと、Pkg もデフォルトでは import されていないので、
import Pkg
と明示する必要があるようだ。

他にも、いろいろとバグがあるようで、例えば、変数のスコープがおかしくなっているようで、
a = 0; for i=1:3; a += i; end
が実行できない。(0.6 では、もちろん実行できる。）

そのうちに、細かいバグが修正されていくのかもしれない。