Error in workspace query in dsytrf for empty matrices

[andreasnoackjensen]

For empty matrices the workspace query returns 0.0 but 0 is an illegal value for lwork. See the following program where I also show the result from dsyev

program test

implicit none
integer*8 :: lwork, info, ipiv(0)
double precision :: a(0,0), w(0), work(1)

lwork = -1_8
call dsyev('N', 'U', 0_8, a, 1_8, w, work, lwork, info)

write(*,*) work

lwork = -1_8
call dsytrf('U', 0_8, a, 1_8, ipiv, work, lwork, info)

write(*,*) work

end program test

[Julien Langou]

Hi Andreas,

Thanks for reporting this. (1) The correct behavior would indeed be to return WORK(1) = 1 when N=0. (2) I have checked a few LAPACK subroutines and unfortunately we are quite inconsistent with respect to this. You reported DSYTRF. DGEQRF has the same issue. DGELS and DSYSV are fine. So I will try to correct DGEQRF and DSYTRF and we will take it from there. Not much point in starting a systematic check.

Cheers,
Julien.

[lawrence mulholland]

This may be down to interpretation.
Routines may say something like "immediate return when N=0"
and may also say "WORK is not referenced N=0".

To be consistent you would have to have a rule as to which takes precedence:
the workspace query or the immediate return.
If the former then you would have to remove the "WORK not referenced comments"
from routines with workspace query mechanisms.

[andreasnoackjensen]

There might be different solutions, but what I wanted to avoid is illustrated in the the program below. If I run a work space query first and then use it I get the error message

** On entry to DSYTRF, parameter number 7 had an illegal value

program test

implicit none
integer :: lwork, info, ipiv(0)
double precision :: a(0,0), w(0), work(1)

lwork = -1
call dsyev('N', 'U', 0, a, 1, w, work, lwork, info)
lwork = work(1)
call dsyev('N', 'U', 0, a, 1, w, work, lwork, info)

lwork = -1
call dsytrf('U', 0, a, 1, ipiv, work, lwork, info)
lwork = work(1)
call dsytrf('U', 0, a, 1, ipiv, work, lwork, info)

end program test