Newton Sketch

SIAM Optimization に「Newton Sketch」というタイトルの論文が出ていたので、チェックをしてみたが、なかなかに面白かった。

普通に Newton 法を実行すると Hessian $$\nabla^2 f(x)$$ が重たいわけだけど、ある構造を満たすような正定値対象行列 $$S$$ を用いて $$S \nabla^2 f(x)^{1/2}$$ のように変更する。
ここで、$$SS$$ は確率的に選んでいて、しかも行列積などが高速に行える構造を持っているようなものを作っている。

このようにすると、stochastic gradient よりも速い収束が得られる、という話だった。

たしかに、stochastic gradient は反復回数が非常に多くて十分な収束が得られないので、こういったように改善する、というのは面白いアイデアだと思う。

論文では、self-concordant なども議論していて、SDP なども解けるとのこと。
ただ、数値実験では SDP が出ていないので、たぶん SDP の構造にはあまりあてはまっていない可能性もある。
そういった可能性もあるが、SDP でどうやって性能を引き出せるか、ということに注目すれば、新しいテーマにできるかもしれない。

みかん最適化分割問題

今朝の NHK で、「みかんの皮を切っていって、酉の形ができあがる」というスゴ技をやっている人がいました。

「あたらしいみかんのむきかた」という本(Amazon の「あたらしいみかんのむきかた」)にはみかんの皮のむき方が書いてあるようなので、実際の絵などは、そちらを見ると分かりやすいかと思います。

これがどうすごいか、というのを数学的な視点からいうと

１．みかんの皮に余りがない（つまり、酉の形になっている皮を戻すと球形に戻る）
２．酉の形は連結である
３．みかんのヘタの部分が酉の目の位置に来ている

の３点が挙げられます。
数理最適化で考えてみると、「３次元多様体から連続的に変化させて得られる２次元多様体で、もっともイメージに近いような変化を求める」というようなことになるかと思われるわけです。

この場合の計算方法としては、次のうちのどちらかなぁ、と思います。
（１）多様体から多様体への写像を求める
（２）皮を離散化して、座標ごとの連結情報を維持しながら、各座標の写像を計算する


うーん、こういったのを頭の中で計算できてしまうとは、世の中には色んな達人がいるものです。

Julia の型変換は、なかなかに難しい (SDPA を Julia で実装する：その３)

Julia では、function の引数に型指定できるので、間違った入力が来ても、引数チェックの時点ではじくことができる。

たとえば、
function func1(x::Int, y::Float64)
とすれば、x は Int, y は Float64 が指定できる。
具体的には、
func1(x::Int, y::Float64) = x+y
等と定義すれば、func1(3, 4.5) で計算できる。


ここで、Int の 3 と Float64 の 4.5 の足し算をすると、Float64 の 7.5 が帰ってくる。
そこで、y のほうに整数を入れても自動的に変換できるだろう、と思うと、実はそうではない。
func1(3,4) とすると、「型が合わない」というエラーが出る。

特に難しいのはベクトルの定義で、例えば
a = [-2; 3; 1]
とすると、Array{Float64} では型が合わないので、基本的には Array{Int}, Array{Float64} のそれぞれで関数を作るか、あるいは float で変換するか、というあたりである。

このあたりを、どうすると簡単に SDPA を使えるようにできるか、検討している。

SDPA を Julia で実装してみる（その２）

SDPA sparse format の読み込みの続きを行った。
Matlab と Julia では、文字列処理が異なるので、単純にコピーするわけにもいかず、このあたりが結構手間がかかる。

Julia では、replace で文字列の置換、split で文字列分解ができて、このあたりは強力な操作ができる。

あと、疎行列については Matlab であればメモリー領域を先に確保しておくことができるが、Julia では同じような機能がない様子。
今は小さいファイルしか読み込んでいないので特に問題ないが、大きい問題になると読み込みに時間がかかりそうだ。

SDPA を Julia で実装してみる（その１）

最近になって、SDPA を Julia で実装しなおしている。
主な理由としては、SDPA のメンテナンスコストを抑えるため。

やはり、SDPAのインストールなどを考えると、コンピュータに関する知識がある程度は必要だし、かといって、OpenBLAS などを毎回コンパイルするのも大変。
そこで、プログラミング言語としても、それなりにスピードがある Julia で実装を行ってみている。


ということで、まずは最初のところとして、SDPA sparse format を Julia に読み込むようにしている。
SDPA では、fscanf の文字処理が強力なので、それに任せているが、Julia だと fscanf を使えないので（libcを経由すれば使えるけど、それだとlibcがないところで使えない）、正規表現で読み込むようにしている。

とりあえずファイルからは読み込めたので、あとは疎行列形式に持たせるようにすれば、OK になるはず。
Matlab に近いインターフェイスにできるだけしたいけど、Matlab だとスピードが出ないところもあるので、そのあたりは適宜修正になりそう。

SOCP に対する多項式時間アルゴリズム

Optimization Online に最近掲載された

http://www.optimization-online.org/DB_FILE/2016/11/5713.pdf

をざっと眺めてみると、SOCP に対する多項式時間アルゴリズムが載っていて、内点法とも近いような、それでいてそうでもないような感じで面白いと思う。

標準形を特定の形式に変形することでアルゴリズムを作るってところは Karmarkar の内点法に似たような感じもあるし、途中で体積の計算をしたりしているところは楕円体法を思い出したりもする。

気になるのは、
（１）実際にプログラムを組んでみたら、どれくらいの性能が出るのか？
（２）内点法よりも数値誤差に強いだろうか？

特に（２）のところは、実行可能内点がない問題がとけるようになれば、ぜひ使ってみたいところ。でも、これはやっぱり難しい？

NASA との共同研究での 3D パッキング問題

INFORMS Annual Meeting で聞いた発表のひとつに、NASA との共同研究で 3D パッキング問題を解いている、というものがあった。

ふつうの問題とは若干違っていて、ひょっとしたら数学的性質がかなり違うのかもしれない。
発表の中では、「厳密解を求めるために、今は整数計画問題として解いている」ということで、定式化されている整数計画問題は比較的オーソドックスに定式化されていた。

これは、なかなかに面白いネタだと思う。
例えば、
１）定式化を工夫して、冗長な情報を除去したり、CPLEX などで高速に解けるようにする
２） Bender's decomposition などが使えるか検討する
３） Lagraingian 緩和, SDP 緩和などの緩和手法で近似値を低コストで得る
などなど、いろんな発展があるんじゃないかと思う。

発表のパワーポイント資料を発表者から頂くことができたので、定式化については検討してみようと思う。