CODEWITHSUNDEEP
fortranlapackblas
Rain
Rain

Reputation: 375

gfortran error: zgesvd in lapack

I was doing an svd decomposition of a square matrix A, with A=U S Vdag, and in the fortran code, the line reads

lwork = -1

call zgesvd( 'A', 'A', A%d, A%d, A%m, A%d, S, U%m, U%d, Vdag%m, Vdag%d,work, lwork, rwork, info )

lwork = int(work(1));deallocate(work); allocate(work(lwork))

call zgesvd( 'A', 'A', A%d, A%d, A%m, A%d, S, U%m, U%d, Vdag%m, Vdag%d,work, lwork, rwork, info )

When I compiled with gfortran, it went through with no error or warning. However when I run the program, it shows error with message:

" ** On entry to ZGESVD parameter number 11 had an illegal value "

I could not figure out what went wrong.

For reference, the definitions of the parameters:

type cmatrix
   integer(4) d
   complex(8), allocatable :: m(:,:)
end type

type (cmatrix) A,U,Vdag

allocate(A%m(dim,dim),U%m(dim,dim),Vdag%m(dim,dim))

A%d = dim; U%m = dim; Vdag%d = dim

real(8) S(dim)

Thanks in advance! Xiaoyu

p.s. It should be mentioned that such a program runs smoothly when compiled with ifort, but gfortran gives an runtime error as shown above

--- Problem solved!

It seems that the issue lies in how ifortran and gfortran allocates memory. I defined in the code USV type which:

type USV
   integer is_alloc
   type (cmatrix) U,V
   real(8), allocatable :: S(:)
end USV

When initializing by 

type(USV) Test_usv(:)

allocate(Test_usv(3)),

the value of is_alloc is 0 using intel fortran compiler, while arbitrary number for gfortran. I need to use this value as a criterion for allocating U V matrices:

if (is_alloc.eq.0) then

   allocate(U%m(dim,dim))

end if

Upvotes: 3

Views: 440

Answers (1)

M. S. B.
M. S. B.

Reputation: 29391

The fundamental problem is not a difference between ifort and gfortran. Your approach to the initialization of the variables isn't valid Fortran. Unless you initialize a variable with a declaration, assignment statement, etc., its value is undefined. One way to fix this would be to add a default initialization to the type definition:

type USV
   integer is_alloc = 0
   type (cmatrix) U,V
   real(8), allocatable :: S(:)
end USV

Another approach would be to not track the allocation status yourself and to rely upon the intrinsic function that Fortran provides for this purpose:

if (.NOT. allocated (U%m) ) allocate(U%m(dim,dim))

P.S. Best practice is not to rely upon specific numeric values for kinds. The kind values are arbitrary and are not necessarily the number of bytes of the type. Some compilers use the number of bytes, others don't. One method to specify the number of bytes and to have portable code is to use types provided by the ISO Fortran environment:

use, intrinsic :: ISO_FORTRAN_ENV

integer(int32) :: d
real(real64), allocatable :: S(:)

The types are named after the number of bits. A list of the available types is in the "Intrinsic Modules" chapter of the gfortran manual.

Upvotes: 1

Related Questions