Help us understand the problem. What is going on with this article?

MKLのLapack95インターフェースを使ってFortranで行列を対角化しよう

Fortran Advent Calendar 2019です。

はじめに

FortranでLapackを使う時によくお世話になるMKLですが、マニュアルを見るとよく「Fortran95インターフェース」というのを見かけると思います。
これを見ると引数がとても減っていて、「なんかよさそう。でもどうやって使うんだろう?」と思う方も多いかと思います。
今回は、このインターフェースをgfortranで使う方法について述べます。
これを使えば、たくさんの引数を使わずに、PythonやJuliaのように簡単に対角化ができます。

以前
「FortranでLapackを使った行列の対角化」
https://qiita.com/cometscome_phys/items/e22d329697e378f6763d
という記事を書きましたが、この時は自前でインターフェースを作りました。
今回はMKLに入っているインターフェースを使います。
gfortranで使う場合には、
https://www.xlsoft.com/jp/products/intel/perflib/mkl/11.2/mkl_userguide_lnx/GUID-C3CFB9D7-EC03-4ECC-8F37-7DAAE2C6C93A.htm
にあるように、まずビルドする必要があります。
ビルドして使えるようになるまで、やってみましょう。

OS

OS: Ubuntu 16.04.6 LTS (Xenial Xerus)
gfortran: gcc version 7.4.0 (Ubuntu 7.4.0-1ubuntu1~16.04~ppa1)

Lapack95インターフェースのインストール

まず、Lapack95インターフェースをgfortran用にビルドします。
MKLは入っているとします。
つまり、

echo $MKLROOT

とした時に、適切にMKLのインストール先が入っているとしましょう。
まず、ホームディレクトリに適当なlapack95のインストール先のディレクトリを作成します。
例えば、

cd ~/
mkdir mkl

などですね。

次に、

cd $MKLROOT/interfaces/lapack95/
make libintel64  INSTALL_DIR=~/mkl/ interface=lp64 FC=gfortran

とすると、いま作ったディレクトリにgfortran用のlapack95インターフェースのライブラリやモジュールがビルドされます。
次に、.bashrcなどにF95ROOTという環境変数を設定します。インストール先のディレクトリの情報です。
例えば、

export F95ROOT=$HOME/mkl

ですね。

これでgfortranでも使えるようになりました。

コンパイル

コンパイル方法は
https://software.intel.com/en-us/articles/intel-mkl-link-line-advisor
でご自分の環境に合わせて調べてみましょう。
Ubuntuであれば、

Use this link line: に

 ${F95ROOT}/lib/intel64/libmkl_lapack95_lp64.a -L${MKLROOT}/lib/intel64 -Wl,--no-as-needed -lmkl_gf_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl

Compiler options:

 -I${F95ROOT}/include/intel64/lp64 -m64 -I${MKLROOT}/include

とあると思います。
これを使います。

lapacktest.f90というコードがあった時、コンパイルは、

gfortran lapacktest.f90  ${F95ROOT}/lib/intel64/libmkl_lapack95_lp64.a -L${MKLROOT}/lib/intel64 -Wl,--no-as-needed -lmkl_gf_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl  -I${F95ROOT}/include/intel64/lp64 -m64 -I${MKLROOT}/include

となります。

サンプルコード

以下のようなサンプルコードを用意しました。

lapacktest.f90
subroutine test()
    use lapack95
    implicit none
    integer::n
    real(8),allocatable::Hs(:,:)
    real(8),allocatable::evec(:)
    real(8),allocatable::Vs(:,:)
    complex(8),allocatable::V(:,:)
    complex(8),allocatable::H(:,:)
    real(8),allocatable::Hr(:,:)
    real(8),allocatable::Hi(:,:)
    complex,parameter::ci = (0d0,1d0)
    n = 2
    allocate(Hs(1:n,1:n))

    allocate(V(1:n,1:n))
    allocate(Vs(1:n,1:n))
    allocate(evec(1:n))

    call random_number(Hs)

    Hs = (Hs + transpose(Hs))/2

    call syev(Hs,evec,'N')
    write(*,*) evec

    call syev(Hs,evec,'V')
    write(*,*) evec

    allocate(H(1:n,1:n))
    allocate(Hi(1:n,1:n))
    allocate(Hr(1:n,1:n))
    call random_number(Hr)
    call random_number(Hi)
    H = Hr + ci*Hi

    H = (H + conjg(transpose(H)))/2

    call heev(H,evec)
    write(*,*) evec

    call heev(H,evec,'V')
    write(*,*) evec



    return
end subroutine

program main
    implicit none
    call test()
end program

以前の記事と同じように、対称行列の固有値と複素数エルミート行列の固有値を求めています。

call syev(Hs,evec,'N')

とすれば、対称行列の固有値が求まり、

call syev(Hs,evec,'V')

とすれば固有ベクトルも求まります。簡単ですね。

lapack95インターフェースを使う時は、

use lapack95

としています。

これで、mklのマニュアルにあるいろんなlapack95インターフェースが使えるようになりました。
もう面倒な引数からはおさらばできますね。

Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away
Comments
Sign up for free and join this conversation.
If you already have a Qiita account
Why do not you register as a user and use Qiita more conveniently?
You need to log in to use this function. Qiita can be used more conveniently after logging in.
You seem to be reading articles frequently this month. Qiita can be used more conveniently after logging in.
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away