[8106] in Athena Bugs
vax 7.3M: recipes
daemon@ATHENA.MIT.EDU (James D. Oliver III)
Wed Sep 4 21:33:28 1991
To: bugs@ATHENA.MIT.EDU
Cc: ellis@ATHENA.MIT.EDU
Date: Wed, 04 Sep 91 21:32:36 EDT
From: James D. Oliver III <oliver@ATHENA.MIT.EDU>
System name: hippocrates
Type and version: CVAXSTAR 7.3M (1 update(s) to same version)
Display type: SM
What were you trying to do?
Use the powell.f subroutines and associated subroutines from the
Numerical Recipes Library.
What's wrong:
I had previously been using my own version of powell.f which I had
tyuped by hand from Numerical Recipes. For convenience in future writing, I
decided to switch over to linking the recipes_f library. A program I had
been successfully using with my version of powell crashed under the library
version, giving me a divide by zero error.
I checked the library source code for powell and found two
differences between my code and the library.
Difference 1 in powell.f
------------------------
My code:
20 continue
call linmin(p,xit,n,fret)
if (abs(fp-fret) .gt. del) then
del=abs(fp-fret)
ibig=i
endif
30 continue
if (2.e0*abs(fp-fret) .le. ftol*(abs(fp)+abs(fret))) return
Library code:
12 CONTINUE
FPTT=FRET
CALL LINMIN(P,XIT,N,FRET)
IF(ABS(FPTT-FRET).GT.DEL)THEN
DEL=ABS(FPTT-FRET)
IBIG=I
ENDIF
13 CONTINUE
IF(2.*ABS(FP-FRET).LE.FTOL*(ABS(FP)+ABS(FRET)))RETURN
Difference 2 in linmin.f
------------------------
My code:
ax=0.e0
xx=1.e0
bx=2.e0
call mnbrak(ax,xx,bx,fa,fx,fb,f1dim)
Library code:
AX=0.
XX=1.
CALL MNBRAK(AX,XX,BX,FA,FX,FB,F1DIM)
I checked with the text of Numerical Recipes and my version of the code is
correct. There may be other errors as I only did a cursory investigation,
but I believe that these are enough to suggest that there may be serious
bugs in the library files.
Please describe any relevant documentation references:
Numerical Recipes, but Press et al.
____________________________
Jim Oliver
oliver@athena.mit.edu joliver@hstbme.mit.edu
oliver%mitwccf.BITNET@MITVMA.MIT.EDU
