calculating inverse matrix using lapack
10,428
From the documentation:
ZGETRI computes the inverse of a matrix using the LU factorization
computed by ZGETRF.
You need to run ZGETRF
first:
program test
Implicit none
integer,parameter::M=300
complex*16,allocatable,dimension(:,:)::A
complex*16,allocatable,dimension(:)::WORK
integer,allocatable,dimension(:)::IPIV
integer i,j,info,error
allocate(A(M,M),WORK(M),IPIV(M),stat=error)
if (error.ne.0)then
print *,"error:not enough memory"
stop
end if
!definition of the test matrix A
do i=1,M
do j=1,M
if(j.eq.i)then
A(i,j)=(1,0)
else
A(i,j)=0
end if
end do
end do
call ZGETRF(M,M,A,M,IPIV,info)
if(info .eq. 0) then
write(*,*)"succeded"
else
write(*,*)"failed"
end if
call ZGETRI(M,A,M,IPIV,WORK,M,info)
if(info .eq. 0) then
write(*,*)"succeded"
else
write(*,*)"failed"
end if
deallocate(A,IPIV,WORK,stat=error)
if (error.ne.0)then
print *,"error:fail to release"
stop
end if
end
Comments
-
someone almost 2 years
I want to calculate inversion of a complex matrix. It occurs to me that lapack contains many routines related to algebraic calculation, so I found the subroutine ZGETRI. Unexpectedly, after compiling the following code with "ifort -o out -heap-arrays test.f90 -mkl" and running the file "out", I got an error
"glibc detected./out:free():invalid pointer: 0x00007fee68f76010***"
followed by memory map and finally "aborted(core dumped)". This is very strange to me and I have no idea where the mistake is. By the way, when some errors emerge not in the compiling process but the running process, is there any way to detect where this error is from?
program test Implicit none integer,parameter::M=300 complex*16,allocatable,dimension(:,:)::A complex*16,allocatable,dimension(:)::WORK integer,allocatable,dimension(:)::IPIV integer i,j,info,error allocate(A(M,M),WORK(M),IPIV(M),stat=error) if (error.ne.0)then print *,"error:not enough memory" stop end if !definition of the test matrix A do i=1,M do j=1,M if(j.eq.i)then A(i,j)=(1,0) else A(i,j)=0 end if end do end do call ZGETRI(M,A,M,IPIV,WORK,M,info) if(info .eq. 0) then write(*,*)"succeded" else write(*,*)"failed" end if deallocate(A,IPIV,WORK,stat=error) if (error.ne.0)then print *,"error:fail to release" stop end if end
-
Ross about 8 yearsThe link to documentation is broken, I suggest netlib.org/lapack/lug/node38.html