g95でLAPACK

g95で実対称行列の固有値を求めるサブルーチンには,
  • dsyev (driver),

  • dsyevd (divide and conquer)

  • dsyevx (expert)

  • dsyevr (relatively robust representations)
といくつかある(LAPACKのユーザガイド及びmanページ参照).以下の簡単な例では,最初の3つの固有値は一致した.dsyevrはe-14〜e-15程度の差があった.

Macには,CBLAS, CLAPACKがAccelerate FrameworkとしてOSに添付されている.g95がインストールしてあれば,

g95 foo.f90 -framework accelerate

でライブラリをリンクできる.とっても簡単.

program eigen
implicit none

character, parameter :: jobz = "V", uplo = "L"
integer, parameter :: n = 3, lda = n, lwork = 3*n-1
integer :: info, i
real*8, dimension(lda,n) :: a
real*8, dimension(n) :: w
real*8, dimension(lwork) :: work

a(:,1) = (/1.0d0, 2.0d0, 3.0d0/)
a(:,2) = (/2.0d0, 3.0d0, -4.0d0/)
a(:,3) = (/3.0d0, -4.0d0, 5.0d0/)

call dsyev(jobz, uplo, n, a, lda, w, work, lwork, info)
print *, "info=", info
do i=1, n
print *, "mode:", i, " eigenvalue=", w(i)
print *, " eigenvector", a(:,i)
end do
end program eigen

コメント

このブログの人気の投稿

静的サイト生成ツールの導入

uptex

TeXのソースを分割してTeXShopでタイプセット