久しぶりに MATLAB を使うことになった。 色々やっていて、ふと 「Bidiagonalization of a matrix based on Lapack interface」 にあるものを試してみようと考えた。
>> a=[1,2;3,4;5,6;7,8]; >> [Q,B,P]=bidiag_lapack (a) 関数または変数 'bidiag_lapack' が認識されません。 |
ああ、MATLAB には bidiag_lapack() があるわけじゃなくて、 それを提供する、ということか。ダウンロードして展開。
% ls -l bidiag_lapack_interface total 240 -rw-r--r--@ 1 mk staff 6054 9 11 2014 bidiag_lapack.m -rw-r--r--@ 1 mk staff 477 9 11 2014 bidiag.m -rw-r--r--@ 1 mk staff 86772 9 11 2014 lapack.c -rw-r--r-- 1 mk staff 1078 4 4 10:59 lapack.log -rw-r--r--@ 1 mk staff 5928 9 11 2014 lapack.m -rw-r--r--@ 1 mk staff 1826 9 11 2014 lbidiag.m -rw-r--r--@ 1 mk staff 1311 9 11 2014 license.txt |
LAPACK のサブルーチンを呼んで、 与えられた行列を二重対角化をするための関数を提供するが、 1つだけ C のプログラム lapack.c というのがあり、 それをコンパイルして使う。
lapack.c の注釈を読む |
/* To build on Windows, run: */ /* mex lapack.c */ /* */ /* To build on any other operating system, run: */ /* mex -ldl lapack.c |
MEX というのは、MATLAB Executable を略して出来た用語で、 “MATLAB上で外部の制御機器を制御するための関数で、 C言語やFortran言語で記述されることが多いです。” だそうです。
cp -pr bidiag_lapack_interface ~/Documents/MATLAB |
>> cd bidiag_lapack_interface >> mex -ldl lapack.c 警告: Xcode is installed, but its license has not been accepted. Run Xcode and accept its license agreement. 次を使用中のエラー: mex サポートされているコンパイラが検出されません。オプションについて は、https://www.mathworks.com/support/compilers を参照してください。 |
あれ?あ、そうか、この Mac は Xcode のインストールを強制されなかったやつだ。 Command Line Tool の cc だけじゃいけないのかな? まあ、Xcode を排除する強い理由もないので (他の Mac にはみなインストールしてある)、 インストールするか。
(しばし Xcode のインストールが終わるのを待つ。)
やり直し |
>> mex -ldl lapack.c 'Xcode with Clang' でビルドしています。 ld: warning: search path '/usr/local/lib' is not a directory ld: warning: search path '/usr/local/lib' is not a directory MEX は正常に完了しました。 |
sudo rm -f /usr/local/lib /usr/local/include sudo mkdir /usr/local/lib /usr/local/include |
>> [Q,B,P]=bidiag_lapack (a) Q = 0.1091 0.8295 -0.3945 -0.3800 0.3273 0.4392 0.2428 0.8007 0.5455 0.0488 0.6979 -0.4614 0.7638 -0.3416 -0.5462 0.0407 B = 9.1652 10.9109 0 0.9759 0 0 0 0 P = 1.0000 + 0.0000i 0.0000 + 0.0000i 0.0000 + 0.0000i 1.0000 + 0.0000i >> Q*B*P'-a ans = 1.0e-14 * 0.0888 0.5329 0 0.0888 0.1776 0.3553 0.0888 0.3553 >> |
まあ、無事に動いたようだ。行列のサイズを大きくして試す。
>> a=rand(100,100); >> [Q,B,P]=bidiag_lapack(a); >> norm(Q*B*P'-a) ans = 2.4587e-14 >> |