2017年6月17日土曜日

TeX のソースコードをテキストに変換

TeX のソースコードをテキストに変換するスクリプトが
https://github.com/subbyte/tex2text/blob/master/tex2text.sh
にあった。

これを使うと、あとで英語の文法チェックのソフトなどにある程度一括でかけることができるので、だいぶ便利になる。

2017年6月16日金曜日

研究ネタになりそうなこと (SIAM Optimization 2017から)

先月参加してきた SIAM Conference on Optimization 2017 で、今後の内容になりそうな内容をメモとして残しておく。

ちなみに、セッション情報などは
http://meetings.siam.org/program.cfm?CONFCODE=OP17
で見ると分かる。
例えば、[MS1-2] はセッション MS1 の中の2番目の発表。

なお、ここにある情報は一部を取り除いてあって、完全版は別途あり。


== 2017/05/22 ==

[MS1-2] The Conic Optimizer in Mosek: A Status Report
Mosek が Ver 8 になって、40% 程度高速になったとのこと。
やはり、問題の再定式化などが要綱。
テスト問題として、以下の論文のものを用いていて、
Reduced-Complexity Semidefinite Relaxation of Optmial Power Flow Problems
http://ieeexplore.ieee.org/document/6692873/
データを作るプログラムがgithubで公開されている。
https://github.com/martinandersen/opfsdr

[MS1-3]
de Klerk の話で、First-order method で最もパフォーマンスが出ないのがどういう
状況なのかを解析するのに SDP を使うっていう話。
基本的には、内積は内積であればどの内積を使ってもよく、重要な情報は Gram 行列に集約される。
思うに、kernel trick のようなものなので、再生ヒルベルト空間の性質が聞いているのでは?
解析では、SDP を解くのではなく、特定のLagrange 乗数を使って KKT 条件を満たしている。
この Lagrange 乗数は、まだ幾何学的な解釈がないようなので、そのあたりも調べたら面白そう。
self-concordant にも似たように拡張していて、現在は内点法について拡張中。
たぶん、inexact Newton だったり、quasi Newton だったりで亜種を作れそう。
論文は Optimization Letter に To appear で、ドラフトとしては、
arXiv:1606.09365

[MS25-2]
SCIP-SDP の話。
SDPA も SDP ソルバーの一つに入れてもらっているけど、
最近の性能向上をさぼっていたためか、ちょっとダメな感じ。
もっと性能を上げるようにしないと。

[MS25-4]
Mixed-Integer Sum of Squares を Julia のパッケージで解く、という話。
Julia のパッケージとしては、MultivariatePolynomial.jl, PolyJuMP.jl,
SumofSquares.jl などがあって、これらを使うと、多項式最適化問題から
SDP 緩和を作るのがとても簡単にできる、ということ。
また、Steering a quadcopter through obstacles では、
Mixed-Integer Sum of Squares で
「障害物を避けながら進む」っていう計算ができていて、
これは計算結果をグラフィック表示で楽しむことができるので、
デモンストレーションなどにいいかもしれない。
(ただ、このデモンストレーションの計算については、論文などにはまだなっていない様子。)

[CP8-4]
サプライチェーンのロバスト版を考えることで SOCP が出てくる、ということ。
また、KKT 条件などから、均衡状態が示されている、ということなので、
SOCP の応用として覚えておくと面白そう。
というか、あっちの応用に組み合わせたら面白い。

== 2017/05/23 ==

[IP-3] Tom Luo の通信に関係する最適化問題の話。
interference の計算に WMMSE という座標降下法の一種を用いていて、
WMMSE はデファクトスタンダードとして用いられている、とのこと。
このときの目的関数はなかなかに複雑になっているので、
あとで細かいところを確認したら面白いかもしれない。

[MS39-1] min c^T x s.t. x in P \cap S で、
P は Ax <= b で S は closed set。
Disjunctive cut を outer-product-free という考えを用いて作られている。
テスト問題には GLOBALLib からとっている。
このGLOBALlib には多項式最適化のベンチマークも含まれるようなので、
いろいろとチェックしたら面白そうだ。

[MS41-2]
arXiv:1702.00526 https://arxiv.org/abs/1702.00526
問題は、min{f(x) | Qx = z, x in Conv(X), z in Z}
Subgradient, Proximal Point Method, Augmented Lagrangian Method,
Bundle Method との関係を調べることで、新しいアルゴリズムを作る。
Proximal SDM (Simplical decomposition method) では、
proximal point method の探索方向の近似を計算していて
column generation と同じ感じで探索領域が広がっていくようだ。
nonlinear SDM-GS (Gauss-Seidel) -ALM の方法を提案していて、
収束などを示しているようだ。
column geration が関係しているから、ほかのところにも使えるかも。

[MS41-3]
内点法ベースで Boyd らが作っているのが Forces というソフトウェア。
Augmented Lagrangian Filter という手法の大域的収束性を示している。
メモリが少ないような組み込み環境でも問題が解ける、ということなので、
に適したアルゴリズムとして考えれば、これは化けるかもしれない。
この話では、空冷システムのようなものの最適化計算を組み込み環境で計算していた。
http://ieeexplore.ieee.org/document/7902123/

[MS40-4]
Monte Carlo ので Adaptive Sampling Method を計算している。

[MS56-1]
Stochastic Mixed Integer Programming の計算を PIPS-SBB という
分子限定法の並列計算で行っている。
大まかに各ノードに割りあてする Coarse-grained と、もっと細かい制御で
並列計算する Fine-grained の2段階割り当てを行う。
ただ、200台ぐらいでパフォーマンスに限度が来る様子。
思うに、何台で最大パフォーマンスになるか、っていう予測があれば役に立ちそう。

[MS56-3]
並列計算の評価基準をどうするべきか、っていう話。
計算時間が遅くても、スケーラビリティーだと高く表示されてしまう危険性がある。

[MS57-3]
Birgin の発表。
ARC (adaptive regularization with cubics)
TRACE (a non-standard trust-region method that uses cubic descent)

f(x+s) <= f(x) - alpha*||s||^3
と最後が3乗になっているのがポイント。
これにより、2乗のときよりも全体の反復回数は抑えられる。
その一方で各反復の計算は2次計画問題を解く必要があって、コストが高い。
思うに、ニュートン法でも反復回数がかかるような問題に有効なのでは?
つまり、first-order < ニュートン法 < cubics のような感じ。

今後の研究では、min f(x) s.t. l <= x <= u のように
上限下限など制約を入れていったときにどうなるか、ということ。

The use of quadratic regularization with a cubic descent condition
for unconstrained optimization
(To appear SIAM J. Optim)
https://www.ime.usp.br/~egbirgin/publications/bmregu.pdf
https://www.ime.usp.br/~egbirgin/sources/bmregu/

[MS57-4]
問題の制約が線形制約な問題に対して、
Gradient Projection Algorithm (GPA) と
Linear Constrained Optimizer (LCO) を繰り返していく
PASA (Polyhedral Active Set Algorithm) というのを作った。
また、 Strong second-order assumption という仮定を導入して、
よりいろんなことを証明している。
Polyhedron への射影は PPPOJ という Hager/Zhang の方法を使う。
また、Dual Active Set Algorithm とも関係があり、
計算の途中で preconditined CG などを使ったりする。
あとは、Grippo のような Non-monotone なども使う。

テスト問題は、Cutest から持ってきている。
結果の比較は IPOPT (デフォルトパラメータ使用)と比較している。
あと、BIGBANK というテスト問題集があり、これには条件数が 10^8 などの
テスト問題が含まれている。
今回の結果は、いかに掲載されている、ということ。
An active set algorithm for nonlinear optimization with polyhedral constraints
https://link.springer.com/article/10.1007/s11425-016-0300-6

=== 2017-05-24===

[IP-4] Renegar の subgradient の話。
Plenary にもかかわらず、概要とか短めに切り上げて、思いっきりコアな数式に突入するあたりが、
うっとりするほど素敵。
通常の subgradient とは違った内容の subgradient なわけだけど、
自己双対錐の性質を使って単位ベクトルに移動してから計算するあたりが面白い。
話の内容としては、Optimization Online に乗っていたのがベース。

[CP12-1]
graph-valued version of total variation っていう話。
全てのラグランジュ乗数が定数になる場合には、実は制約なし最適化問題を解けば十分、
というのがとても面白い性質だった。
ちなみに、出てきている問題は2次計画問題に帰着していたけど、これが凸ではないので、
そう簡単には解けないはず。

[CP15-5]
カラテオドリの一般化の話。
このあたりは、理論的な枠組みと実際の数値計算がずれるかもしれないので、
もしずれるのであれば、それをどうやって修正するか、とか考えたら実用に使えるかも。

[CP15-6]
Conic Programming の実行可能でないとかいろいろと分岐するのを
ADMM  の計算過程に基づいて判別している。
このあたりは、SDPA の判別基準とも似たようなことになっているので、
その二つの論文を突き合わせて比較すると、もう少し詳細に踏み込めるはず。
(たぶん、SDPA のパラメータの omega についての考察とか。)
そういった情報をうまく取り入れていくと、面白いことがわかるかも。



2017年5月16日火曜日

IguanaTex の新しいバージョンと新しいバグ

IguanaTex のバージョンが新しくなっていて、なんとベクター画像を扱えるようになっていた。
いままでは基本的に png で画像ができて、それを貼っていたわけだけど、今度はベクターになるので、拡大や縮小したときにもスムーズにできる。


で、そのベクター画像の機能だけど、変なバグがあるようで、カラーにするとおかしなことが起こる。
たとえば、
\textcolor{red}{ここは赤色}
とかしても、なぜかカラーにならない。
ところが、Vector ではなくて Bitmap にするとちゃんと赤い文字になる。
そのあとでVectorとして画像を生成させると、今度はベクター画像で赤い文字になる。
なかなかに面白い現象だ。

いずれにしても、IguanaTex は beamer の100倍程度の便利さなので、いつまでもなくならないでほしいところだ。

Stochastic subgradient descent method for large-scale robust chance-constrained support vector machines

新しい論文
Stochastic subgradient descent method for large-scale robust chance-constrained support vector machines
https://link.springer.com/article/10.1007/s11590-016-1026-4
について読んでみた。

使っている手法としては、よくある感じのものだった。
ただ、理論的な解析がほとんどないので、どうしたんだろうと思っていたら、解いている最適化問題が非凸最適化問題だった。
したがって、first-order method によくあるような O(1/epsilon) などの収束レートは出てこない。

数値実験を見てみると、非凸最適化の定式化でも、そこそこ最適値を得られるようで、局所的最適解にはちゃんと到達できているのかもしれない。
そのあたりを詳しく解析してみると面白そう。

2017年5月2日火曜日

A doubly inexact interior proximal bundle method for convex optimization

Optimization Online に面白そうな内容が載っていたので、チェック。

A doubly inexact interior proximal bundle method for convex optimization
http://www.optimization-online.org/DB_HTML/2017/04/5960.html

基本的なアイデアとしては、inexact  Level-Set method と内点法を組み合わせたところ。
この中では、最適解に収束する、ということが証明されているが、こういった inexact の解法についてはどれくらいの収束レートになるか、ということを扱っている本があるので、それと合わせると収束レートが分かって面白いかもしれない。
気になるのは、内点法のように超一時収束できるのか、あるいはinexact が効いてあまり収束レートは良くないのか、というところ。

あと、主問題の解を構成していたり、Bregman distance を使っていたりするところは、他にも使ってみたい。

2017年5月1日月曜日

Model Building in Mathematical Programming

Model Building in Mathematical Programming
http://as.wiley.com/WileyCDA/WileyTitle/productCd-1118443330.html
をぱらぱらとチェックしてみた。

タイトルにある通り、数理最適化問題にどのように定式化するか、というのがメインだけど、いろんな定式化が載っていて為になる。

たとえば、x = {0} \cup [a, b] \cup {c} のように x の範囲がとびとびになってしまう場合に、どうやると整数計画問題にできるか、など書いてあったりもする。

実際の問題を数理最適化問題にするときには、いろんな定式化があって、どの定式化を選ぶかで計算時間に大きな影響があるので、こういうのを読んでおくっていうのは、いいことなんだろうなぁ、と思ったりする。

というか、日本語訳を出したらいいんじゃないだろうか、とも思う。
できれば最新の第5版で。(いろいろと追加されているし。)


2017年4月26日水曜日

Julia の Windows 用エディタとしては Notepad++ がベスト

Julia の Windows 用エディタとして、Atom や Visual Studio Code, IJupyter などを試してきたが、これらでは今一つであった。
最大の問題点は、動作が遅いことで、特に Atom はショートカットをクリックしてから30秒以上しないと Julia のソースファイル編集に入れない状態だった。

結局のところ、Notepad++ と REPL を使うのが現状ではベストのようだ。

便利なところとしては、

1.Julia 用の表示にできる

https://github.com/JuliaLang/julia/blob/master/contrib/Julia_Notepad%2B%2B.xml
このファイルをダウンロードして、Notepad++ のフォルダに入れて、[言語]からユーザー定義として「インポート」すれば、Julia 用にできる。

2.編集中のファイルのあるフォルダで Julia の REPL を起動できる。
実行コマンドの登録は、
http://so-zou.jp/software/tech/tool/editor/notepad-plus/commands/
に載っている方法でできるが、Julia の REPL を起動するときには、

cmd /c "cd /d $(CURRENT_DIRECTORY) & c:\Users\ユーザー名\AppData\Local\Julia-0.5.0\bin\julia.exe"

を登録すればよい。ただし、ユーザー名と Julia のバージョンで若干の修正が必要である。
このショートカットを、例えば Alt+F5 に割り当てれば、Alt+F5 で REPL が起動する。
REPL の中で cd で移動するのは大変なので、これが地味に便利である。

3.自動保存ができる。
Notepad++ のプラグインで Autosave を入れると、オプションで設定すれば Notepad++ からフォーカスが外れたタイミングで自動保存をできる。
つまり、Notepad++ で example1.jl を編集しているときに、Alt+Tab で REPL のウィンドウに移動したときに自動的に example1.jl が保存される。
そこで、include("example1.jl"); を実行すればよい。
(何回も実行するのであれば、ctrl+p で前のコマンドを REPL は呼び出せるので、Alt+Tab, ctrl+p, Enter で実行できる。)

4.速い。
Notepad++ を起動してから Julia のファイルの編集に入るまでに1秒以内でできる。自分のところでは、Atom の30倍以上高速である。


とりあえず、Notepad++ を使い込んでみようと思う。



2017年4月25日火曜日

scip には、合わないデータだったみたい

昨日までだいぶかかってインストールした SCIP に早速データを投入して数値実験。
結果としては、CPLEX のほうが性能が良かった。
このあたりは入力データ依存なので、たまたま SCIP に合わないデータなのかもしれない。

分かった点としては、
1.primal bound (実行可能点)が5時間ぐらいたっても1つも得られない。したがって duality gap も出てきていない。
   CPLEX の場合は、数分で実行可能点をとりあえずは発見している。ついでに、dual bound も CPLEX の値のほうがいい。
2.メモリ使用量が劇的に増える。
   分枝の方法があってないのかもしれないけど、とりあえず5時間ぐらいで300万ノードぐらいが生成されていて、100GB のメモリを使用。

このままだと「コンピュータに負担はかかるけど、有益な情報が得られない」ということになりそうなので、計算をストップすることにしてみた。

やはり、入力データにどれだけ性能が出せるかはやってみないと分からないものだ。

2017年4月24日月曜日

SCIP を Linux 上の Matlab から呼び出す、ってお話

SCIP を Linux 上 の Matlab から呼び出すのは、あまり推奨されていないのかもしれないが、コンパイルを頑張って、とりあえずはサンプルプログラムまでは実行できるようになった。(ちなみに、Matlab から呼び出すのでなく単純にコマンドラインから使うだけなら、それほど険しい道のりではない。)

ここにはインストールの記録を残しておく。なお、Ubuntu 系の Linux を想定しているので、ディレクトリ名などは適宜変更が必要かもしれないし、apt-get で必要なファイルなどをインストールする必要もあるかもしれない。
また、インストールするディレクトリは、とりあえず $HOME/lib としている。

1. SoPlex をコンパイル
cd $HOME/lib
wget http://soplex.zib.de/download.php?fname=soplex-3.0.0.tgz
tar xzf soplex-3.0.0.tgz
cd soplex-3.0.0
make SHARED=true

(SHARED=true は Matlab の mex にリンクするときに必須なので忘れないように。)

2. ZIMPL をコンパイル
cd $HOME/lib
wget http://zimpl.zib.de/download/zimpl-3.3.4.tgz
tar xzf zimpl-3.3.4.tgz
cd zimpl-3.4.4
make SHARED=true

(このあたりの wget は、ライセンスの関係でうまくいかないかもしれないので、ホームページからダウンロードしたほうがいいかもしれない。)

3. IPOPT をコンパイル
cd $HOME/lib
wget  https://www.coin-or.org/download/source/Ipopt/Ipopt-3.12.7.tgz
tar xzf Ipopt-3.12.7.tgz
cd Ipopt-3.12.7
./configure
make

(make install は特に必要なし。)

4. SCIP をコンパイル
cd $HOME/lib
wget http://scip.zib.de/download.php?fname=scip-4.0.0.tgz
tar xzf scip-4.0.0.tgz
cd scip-4.0.0
make SHARED=true

(このあたりのコンパイルで、ZIMPLE や SOPLEX, IPOPT がどこにあるか指定する必要あり。画面上のメッセージは分かりにくいけど、20回ぐらいトライ&エラーをすれば、なんとかなった。)

5. OPTI をコンパイル
cd $HOME/lib
wget https://github.com/jonathancurrie/OPTI/archive/master.zip
unzip master.zip
cd OPTI-master/Solvers/Source/scip

ここで、以下のようなスクリプトでコンパイル。

******************************
#!/bin/sh -x

SCIP_DIR=$HOME/lib/scip-4.0.0
SOPLEX_DIR=$HOME/lib/soplex-3.0.0
OPTI_DIR=$HOME/lib/OPTI-master
IPOPT_DIR=$HOME/lib/Ipopt-3.12.7
mex -largeArrayDims -output scip.mexa64 -I$SCIP_DIR/src -I$SOPLEX_DIR/src -I$OPTI_DIR/Solvers/Source/scip/Include \
    -I$IPOPT_DIR/Ipopt/src/Common -I$OPTI_DIR/Solvers/Source/opti\
    -L$SCIP_DIR/lib/shared \
    scipmex.cpp scipeventmex.cpp scipnlmex.cpp -lscip -lnlpi.cppad -lut
******************************

基本的には、リンクで見つからない関数があれば、そのつど nm コマンドで *.so を全解列挙して、grep で検索すればよい。

6. Matlab でパスを追加
addpath(genpath(getenv('HOME'),'/lib/OPTI-master'));

7. Matlab を以下のスクリプトで起動

******************************
export LD_LIBRARY_PATH=/usr/lib/gcc/x86_64-linux-gnu/6.3.0:$HOME/lib/
scip-4.0.0/lib/shared:/usr/lib/x86_64-linux-gnu:$HOME/lib/scip-4.0.0/
lib/shared:$LD_LIBRARY_PATH

LD_PRELOAD=/usr/lib/gcc/x86_64-linux-gnu/6.3.0/libgcc_s.so:/usr/lib/gcc/x86_64-linux-gnu/6.3.0/libstdc++.so:/usr/lib/gcc/x86_64-linux-gnu/6.3.0/libgfortran.so:$HOME/lib/scip-4.0.0/lib/shared/libscip.so:/usr/lib/x86_64-linux-gnu/libreadline.so:$HOME/lib/scip-4.0.0/lib/shared/libnlpi.cppad.so:$HOME/lib/scip-4.0.0/lib/shared/libscipsolver.so  /usr/local/bin/matlab -nodisplay $*
******************************

実行時にうまく行かない場合も、nm と grep で探し当てて、LD_LIBRARY_PATH と LD_PRELOAD に載せれば何とかなる。
(なぜか LD_PRELOAD だけでは、うまく行かない時があるので、LD_LIBRARY_PATH にもフォルダをあわせて書く方がよいようだ。)

7. サンプルソースの実行
あとは、$HOME/lib/OPTI-master/Examples/Global_NonlinearProgramming.m にある SCIP 用のサンプルを実行すればよい。
具体的には、Example 4 - White Box Quartic となっているところである。



途中途中であまりスリムになっていないところもあるが、これでサンプルプログラムを実行できる。
Matlab から $ help scip とすれば使い方も分かるようになっている。

たぶん、Mac でも同じような手順でコンパイルできるかと思うが、Mac は OS 自体が複雑なので、もっと大変なのかもしれない。


2017年4月20日木曜日

SCIP がインストールできずに、てこずる話

どうやら SCIP というのが性能がいいらしい、という話を今更ながらに聞いてきたので、今日インストールしようとし始めた。

扱いたいデータは Linux サーバーのMatlab 上にあるので、SCIP を Linux の Matlab から呼び出すようにしないといけない。
どうやら、YALMIP を使って OPTI 経由で SCIP を呼び出すのが流儀のようなので、インストールしてみる。

今日のところインストールしたのは、
SCIP
Soplex
ZIMPL
YALMIP
OPTI
Ipopt
まで。

で、OPTI は Windows 版はあるけど Linux 版がないので、scip のインターフェイスを直接 mex でコンパイルしようと模索。
すると、なんと libscip.a が -fPIC がついていなかった。

明日のうちに SCIP のコンパイルからやり直さないと。。。

できれば、SCIP が apt-get でインストールできると便利ではあるけど、apt-get でインストールできるようになるのは、現状では難しいと思う。

2017年4月12日水曜日

救急車再配置の実例

今朝のニュースで、リアルタイムでの救急車再配置を実際に行っている例が取り上げられていた。
前に調べた時には、日本では実例が見つからなかったので、大変参考になった。
実際にどのように行っているのか、情報がどこかにあるか調べてみようと思う。


2017年4月7日金曜日

IT Text 数理最適化

世の中には、非常に不思議だと思うこともあって、「IT Text 数理最適化」の Amazon のページにレビューが1件もなかったりする。

数理最適化の基礎を勉強するうえでとてもいい本だと思うし、自分もちょくちょく参考にしていたりするし、他の人から「数理最適化って、どの本で勉強したらいいですか?」って聞かれたら「この本がおススメ」って伝えていたりもするし。

なんでレビューが1件もないのか謎だ。不思議だ。



2017年3月31日金曜日

TeX で使われていない数式番号を見つけ出す

TeX で \label をつけて数式番号を振ってあるけど、その数式番号をどこでも使ってなかった、ということが、よくある。
とりあえずは、\label がついているものすべてに \ref が使用されているかどうかを見ればわかるので、Linux なら次のようなスクリプトで \label と \ref の diff を取ることができる。

(ただし、paper.tex は適宜変更のこと)
#!/bin/sh
grep -o -e "ref{[^}]*}" paper.tex | sort | uniq  | sed "s+ref{++" | sed "s+}++" > __ref.txt
grep -o -e "label{[^}]*}" paper.tex | sort | uniq  | sed "s+label{++" | sed "s+}++" > __label.txt

diff __ref.txt __label.txt
rm __ref.txt __label.txt

2017年3月3日金曜日

Optimality condition and complexity analysis for linearly-constrained optimization without differentiability on the boundary

Optimization Online で

Optimality condition and complexity analysis for linearly-constrained optimization without differentiability on the boundary
Gabriel Haeser, Hongcheng Liu, Yinyu Ye
http://www.optimization-online.org/DB_HTML/2017/02/5861.html

というのが掲載されていたので、チェックをしてみた。
基本的には、KKT 条件を epsilon 誤差だけ許容するような条件に拡張した場合にどうなるか、という点が議論されていた。

実際の数値計算では、厳密に(誤差0で)KKT条件を満たせることは、ほぼないので、epsilon-KKT 条件のようなものが計算されているはずで、その理論展開が普通の KKT 条件と似たようなところで展開されている、というのが面白い。
つまり、「KKT 条件が成り立たないはずの数理最適化問題が、KKT条件を想定して作られたソフトウェアでも、なぜか答えが求まっちゃう」ということが実際には良くあるわけだが(例えば、Slater 条件を満たさない SDP の解を SDPA が求めちゃったり)、こういう研究を見てみると、そういう現象の背景が推察できたりする。

2017年2月28日火曜日

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 でどうやって性能を引き出せるか、ということに注目すれば、新しいテーマにできるかもしれない。



2017年1月16日月曜日

みかん最適化分割問題

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

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

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

1.みかんの皮に余りがない(つまり、酉の形になっている皮を戻すと球形に戻る)
2.酉の形は連結である
3.みかんのヘタの部分が酉の目の位置に来ている

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

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


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


2017年1月13日金曜日

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

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 を使えるようにできるか、検討している。




2017年1月11日水曜日

SDPA を Julia で実装してみる(その2)

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

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

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

2017年1月10日火曜日

SDPA を Julia で実装してみる(その1)

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

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


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

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