Why i get zero results in LAPACK

Help
1

Fortran 90 CODE

Implicit None

Integer::n,nrhs,lda,ldb,info

Real,Allocatable::a(:,:),b(:),x(:)

Real,Allocatable::ipiv(:)

n = 3

Allocate(A(n,n),b(N),x(N),ipiv(n))

A(1,1)=3

A(1,2)=2

A(1,3)=-1

A(2,1)=2

A(2,2)=-2

A(2,3)=4

A(3,1)=-1

A(3,2)=1/2

A(3,3)=-1

b(1) = 1

b(2) = -2

b(3) = 0

x = b

nrhs = 1

lda = n

ldb = n

call dgesv(n,nrhs,A,lda,ipiv,x,ldb,info)

Write(*,*)x(1),x(2),x(3)
End

gfortran test.f90 -o test liblapack.a libblas.a libtmglib.a

0.00000000 0.00000000 0.00000000

2

You have to use double precision real numbers with LAPACK routines which begin with ‘d’.

Try calling sgesv instead.

3

thank you for your replying
I used double precision but I got Incorrect Results
-0.50000000000000033 1.5000000000000007 0.50000000000000033

4

This is not what you think it is. Fortran uses integer division (similar to how Python 2 did), so this evaluates to 0. Use 0.5 or 1.0/2 instead.

1 Like
5

thank you very much Mr, I used 0.5 and the problem was solved.

6

Note a few things:

  • 0.5 is a single precision constant, and when you assign it to a double precision variable you may not have all the precision you want. If you’re lucky, like with 0.5, you have the full precision, because 0.5 has an exact representation in both single and double precision. But with 0.3 for instance you will be unlucky because it does not have an exact representation. So, when you want to make calculation in double precision, take the habit of writing double precision constants: 0.5 shoud be 0.5d0.
  • ipiv in you code is supposed to be an integer array.
  • you can take advantage of the array syntax to simplify a lot of things:
program test
implicit none

integer :: n,nrhs,lda,ldb,info
double precision, allocatable :: a(:,:),b(:),x(:)
integer, allocatable :: ipiv(:)

n = 3
allocate(A(n,n),b(N),x(N),ipiv(n))

A(1,:) = [  3d0,   2d0, -1d0 ]
A(2,:) = [  2d0,  -2d0,  4d0 ]
A(3,:) = [ -1d0, 0.5d0, -1d0 ]

b(:) = [ 1d0, -2d0, 0d0 ]

x = b

nrhs = 1
lda = n
ldb = n

call dgesv(n,nrhs,A,lda,ipiv,x,ldb,info)

write(*,*) x(:)
end
1 Like