Question

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 
Was it helpful?

Solution

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
Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top