本篇博文简要介绍使用MKL函数库计算方阵的逆矩阵
代码如下:
program MKL_getrfANDgetri
use lapack95
implicit none
integer, parameter :: n = 3
integer :: i, j, ipiv(n)
real(kind=8) :: a(n,n), aa(n,n)
call random_seed()
call random_number(a)
aa = a
write(*,'(1x,a)') "a = "
do i = 1, n
write(*,'(*(f12.6,3x))') a(i,:)
end do
!// 使用lapack求解逆矩阵
call getrf( a, ipiv )
call getri( a, ipiv )
write(*,'(1x,a)') 'inv(a) = '
do i = 1, n
write(*,'(*(f12.6,3x))') a(i,:)
end do
write(*,'(1x,a)') "checking..."
aa = matmul(aa,a) !// 原矩阵与其逆矩阵的结果为单位矩阵
do i = 1, n
write(*,'(*(f12.6,3x))') aa(i,:)
end do
end program MKL_getrfANDgetri
结果如下:
a =
0.959299 0.268247 0.274620
0.013673 0.082084 0.275984
0.056097 0.730892 0.177709
inv(a) =
1.072178 -0.876915 -0.295021
-0.074784 -0.888503 1.495423
-0.030876 3.931107 -0.430155
checking...
1.000000 0.000000 -0.000000
-0.000000 1.000000 0.000000
-0.000000 0.000000 1.000000
这里值得注意的是,在调用getri计算逆矩阵之前,应该先调用getrf先进行矩阵的LU分解