LAPACKの利用 #

1. 実行環境 #

※WSLはWindows Subsystem for Linux、PSはPowershellを意味します。

2. LAPACKのインストール #

英語では次の説明が参考になります。

LAPACK for Windows

日本語では下記が参考になります。

WindowsのMinGWでlapackを使う

ただし、日本語の上記説明だけではインストールができず、以下を変更する必要がありました。

• 日本語説明の11.を変更 Specify native compilers にて、LAPACK for Windowsに従い下記のように変更しました。

  • Compilers
    • C ………. C:/TDM-GCC-64/bin/x86_64-w64-mingw32-gcc.exe
    • C++ …….. C:/TDM-GCC-64/bin/x86_64-w64-mingw32-g++.exe
    • Fortran …. C:/TDM-GCC-64/bin/x86_64-w64-mingw32-gfortran.exe

• 日本語説明の12.を変更 LAPACK for Windowsに従い下記のように変更しました。

  • BUILD_SHARED_LIBSのオプションをON
  • CMAKE_GNUtoMSのオプションをON

ここまででLapackが使用できるようになりました。

3. LAPACKを含むFortranコードをPythonで実行する #

3.1. 準備 #

下記のファイル構造をしているとします。

.
├── inverse_matrix.f90
└── main.py

下記リンクより、上記のファイルを入手できます。
lapack_test.zip (約1kB)

subroutine gen_inv(N, A)
    implicit none
    integer, intent(in) :: N
    double precision, intent(inout) :: A(N, N)

    integer :: lda, ipiv(N), info
    integer :: lwork
    double precision :: work(N)

    lda = N
    lwork = N
    call dgetrf(N, N, A, lda, ipiv, info)
    call dgetri(N, A, lda, ipiv, work, lwork, info)

    return
end subroutine gen_inv

import numpy as np
import os

os.add_dll_directory("C:/TDM-GCC-64/bin")
import invmat_f90  # type: ignore


if __name__ == "__main__":
    print(invmat_f90.gen_inv.__doc__)

    A = np.zeros((2, 2), dtype=np.float64)
    for i in range(2):
        for j in range(2):
            A[i, j] = i + j + 4.3

    print(A)
    A_inv = invmat_f90.gen_inv(a=A, n=A.shape[0])
    print(A_inv)

    print(np.dot(A, A_inv))

3.2. Fortranをコンパイル #

f2py -c inverse_matrix.f90 -m invmat_f90 -lblas -llapack

上記を実行し、成功すると環境にもよりますが次のpydファイルが得られます。

.
├── inverse_matrix.f90
├── invmat_f90.cp313-win_amd64.pyd
└── main.py

3.3. Pythonの実行 #

実行すると、次の結果を得ます。

PS C:\work\f2py\test\lapack_test> python .\main.py
a_inv = gen_inv(a,[n])

Wrapper for ``gen_inv``.

Parameters
----------
a : input rank-2 array('d') with bounds (n,n)

Other Parameters
----------------
n : input int, optional
    Default: shape(a, 0)

Returns
-------
a_inv : rank-2 array('d') with bounds (n,n)

[[4.3 5.3]
 [5.3 6.3]]
[[-6.3  5.3]
 [ 5.3 -4.3]]
[[1.00000000e+00 5.32907052e-17]
 [1.63424829e-15 1.00000000e+00]]
PS C:\work\f2py\test\lapack_test>

始めの print(invmat_f90.gen_inv.__doc__) によって得られた結果を確認すると、f2pyによって自動的に生成された関数の引数を確認することができます。

Fortranではsubroutine gen_inv(N, A, A_inv)と引数が3つありましたが、f2pyによって実際にpythonから呼び出す場合はAのみで良くなっており、返り値としてA_invが登録されていることが分かります。Aのサイズは自動的にoptionと判断されており、指定しなければA.shepe[0]から決まることが書かれています。

つまり、最小の引数の数で呼ぼうとすれば、A_inv = invmat_f90.gen_inv(A) だけで良いということです。

4. tips #

少し触っていましたが、f2pyを利用したFortranの実行は癖がいくつかあります。できるだけ手堅い動作をさせるためには、次の注意に従うとよいです。

• Fortranを書く時の注意

  • intent(in), intent(out)のみを利用して入出力を管理すること。intent(inout)の挙動がわかりにくい。
  • 引数の順番の基準がわかりにくいため、関数、サブルーチンの引数は全て名前の指定して渡すこと。
    • ret = f90_func(a=a_py, b=b_py, c=x_py)と全ての引数に=をつけて渡すのが良いでしょう。

• Pythonを書く時の注意

  • os.add_dll_directory("C:/path/to/your/librarydll")を追加すること。
    ただし生成されたpyd内でlibgfortran-5.dll や libquadmath-0.dllに依存する関数をしていない場合は不要。

    次のコードのように、f2pyで生成したpydをimportする前に書いておかなければ、ImportErrorが出ます。

    もしpyd内でlibgfortran-5.dll や libquadmath-0.dllを利用していると、f2pyによって生成されたpydを

    Python

    import pyd_name_generated_by_f2py

    とimportすると、次のエラーが出ます。

    Txt

    ImportError: DLL load failed while importing xabc: 指定されたモジュールが見つかりません。

    これを回避するためには

    Python

    import os
    os.add_dll_directory("C:/TDM-GCC-64/bin")
    import pyd_name_generated_by_f2py

    として、明示的に指定する必要があります（私の環境では、C:/TDM-GCC-64/bin 以下にlibgfortran-5.dll や libquadmath-0.dllが存在しているためです）。

    セキュリティの問題から、ユーザー環境変数PATHまたはシステム環境変数PATHに “C:/TDM-GCC-64/bin” が書かれていても参照されなくなったようです。