H kan lane 2008-01-26 11:52:37
We use LAPACK in a Finite Element program for structural mechanics
applications. It has been and still is a task of ours to port it from
Unix to a Linux cluster. However, we’ve stumbled upon some real
The code currently gets infinitely “stuck” in the dlarfg routine. To be
even more precise, this code is similar to where it enters and never
BETA = -SIGN( DLAPY2( ALPHA, XNORM ), ALPHA )
SAFMIN = DLAMCH( ‘S’ ) / DLAMCH( ‘E’ )
IF ( ABS( BETA ).LT.SAFMIN) THEN
DO WHILE ( ABS( BETA ).LT.SAFMIN )
RSAFMN = ONE / SAFMIN
KNT = 0
KNT = KNT + 1
CALL DSCAL( N-1, RSAFMN, X, INCX )
BETA = BETA*RSAFMN
ALPHA = ALPHA*RSAFMN
Beta is forever equal to a small number and SAFMIN to 1. Have you
encountered this problem before and can you offer any advice? Thanks in
H kan Lane
(Remove CUT., SPAM. and PLEASE. for email address)
H kan lane 2008-01-26 11:52:40
Now the email address is correct in the header.
Alois steindl 2008-01-26 16:20:40
H kan Lane
I guess that your version of dlamch isn’t working properly.
On an i686 Linux system I get
Epsilon = 2.22044605E-16
Safe minimum = 2.22507386E-308
so the ratio safmin should be _very_ small (~1.0e-294), not close to
For safmin=1 the above loop will never end.
Which Lapack version are you using?
H kan lane 2008-01-26 16:20:49
Thank you for your reply!!
The header comment says that it was constructed in March 26, 1990.
Is that the same as the one you’ve used?
Alois steindl 2008-01-26 16:20:50
H kan Lane
1990 sounds quite old for a Lapack release.
If you are working with Linux, you should be able to find precompiled
versions. (Precompiled and optimized versions exist also for other
Unix vendors, like IBM AIX, HP, …) I would also suggest to look at
the package atlas, it is drastically faster than plain Lapack/Blas.
Intel also provides optimized Lapack packages for Windows and Linux.
I would suspect that the person, who build your version, didn’t know
Fortran extremely well and therefore didn’t get the prebuilding step
right. Maybe he didn’t know, what dlamch is intended to return.
H kan lane 2008-01-26 16:20:52
That’s funny. If I split the two routines into one each, I get
eps=3.4916190425559366E-003 and the safe minimum =
1.0404122252830053E-268. Hence, SAFMIN becomes a very high number. The
version I’m using has been downloaded from the web and seems to be the
How would I go about replacing it?
Thank you again.
H kan lane 2008-01-26 16:20:54
Sorry, it should be the other way around.
Also, dlamch won’t work at all unless I insert some print statements.
Excuse the confusion.
Alois steindl 2008-01-26 16:20:55
H kan Lane
you are not giving any details, so its quite hard to give any
meaningful answer. From your first posting I conclude, that you are
Now you would have to tell, which machine you are running on, which
Linux (or other OS) version, which Fortran compiler.
Also “From the web” is _very_ unspecific. The best place to look is at
your vendor’s (eg. Linux or Fortran’s distributor) homepage.
Depends, how it is installed. On a Redhat system one would just say
It is also quite simple – but maybe not that much efficient – to
compile it yourself using the source files on www.netlib.org.
H kan lane 2008-01-26 16:20:57
The compiler version is pgf90. It runs on a cluster of Redhat Linux 8.0.
This is installed on PC:s.
It’s apparent that dlamch is the “culprit”. I’ve tried to download it to
no avail. It still won’t sork correctly.
Thanks once more.
It’s now apparent tht
Alois steindl 2008-01-26 20:49:08
H kan Lane
now we are certainly quite off topic for this group, so I think thats
my last posting in this thread.
If you have a close look at the Redhat cd’s, you will find lapack and
We know, that dlamch behaves strangely on your system, but we don’t know
why. Possible reasons are: (certainly not complete)
I have no idea, how good pgf90 is, but do you really need a
f90-compiler? For Lapack and Co g77 is fine. You could also try the
free Intel compiler or the comercial products from Lahey or Nag.
Dlamch comes with Lapack, where did you try to download it?
H kan lane 2008-01-26 20:49:20
Please reply to my email address without CUT., SPAM. and PLEASE. if you
wish. Thank you greatly.
I don’t know whether we have those compilers. Are they freely available
over the web?
I believe I downloaded it from the www.NetLib.org site.