This bug occurs on gfortran 4.1 and 4.2 . I think it is not a gfortran specific bug; I checked g77 and g95 on gcc 3.4.6.. I had compiled a legacy fortran77 code and foud a bug; $ gfortran -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali Aligned length= 89, RMSD= 6.41, TM-score= 0.24257, ID=0.042 $ gfortran -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali Aligned length= 91, RMSD= 6.35, TM-score=0.24762 , ID=0.024 $ pgf77 -O1 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali Aligned length= 91, RMSD= 6.35, TM-score=0.24762, ID= 0.024 $ pgf77 -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali Aligned length= 91, RMSD= 6.35, TM-score=0.24762, ID=0.024 That optimization error is on a bellow subroutine; ******************************************************************** * Dynamic programming for alignment. * Input: score(i,j), and gap_open * Output: invmap(j) ******************************************************************** SUBROUTINE DP(NSEQ1,NSEQ2) PARAMETER(nmax=5000) LOGICAL*1 DIR common/dpc/score(nmax,nmax),gap_open,invmap(nmax) dimension DIR(0:nmax,0:nmax),VAL(0:nmax,0:nmax) REAL H,V *** initialize the matrix: val(0,0)=0 do i=1,nseq1 dir(i,0)=.false. val(i,0)=0 enddo do j=1,nseq2 dir(0,j)=.false. val(0,j)=0 invmap(j)=-1 enddo *** decide matrix and path: DO j=1,NSEQ2 DO i=1,NSEQ1 D=VAL(i-1,j-1)+SCORE(i,j) H=VAL(i-1,j) if(DIR(i-1,j))H=H+GAP_OPEN V=VAL(i,j-1) if(DIR(i,j-1))V=V+GAP_OPEN IF(( D.GE.H).AND.(D.GE.V)) THEN DIR(I,J)=.true. VAL(i,j)=D ELSE DIR(I,J)=.false. if(V.GE.H)then val(i,j)=v else val(i,j)=h end if ENDIF ENDDO ENDDO *** extract the alignment: i=NSEQ1 j=NSEQ2 DO WHILE((i.GT.0).AND.(j.GT.0)) IF(DIR(i,j))THEN invmap(j)=i i=i-1 j=j-1 ELSE H=VAL(i-1,j) if(DIR(i-1,j))H=H+GAP_OPEN V=VAL(i,j-1) if(DIR(i,j-1))V=V+GAP_OPEN IF( V.GE.H) THEN j=j-1 ELSE i=i-1 ENDIF ENDIF ENDDO c^^^^^^^^^^^^^^^Dynamical programming done ^^^^^^^^^^^^^^^^^^^ return END Files; http://zhang.bioinformatics.ku.edu/TM-align/TMalign.f http://zhang.bioinformatics.ku.edu/TM-align/benchmark/1aquA.pdb http://zhang.bioinformatics.ku.edu/TM-align/benchmark/1avaC.pdb
Created attachment 13727 [details] A legacy fortran77 program
Created attachment 13728 [details] A input data file of TMalign
Created attachment 13729 [details] A input data file of TMalign
Thank you, Tobias I had missunderstood the default optimization level for gfortran but the issue exists, I think. I had traced side effects of optimization levels for the legacy program; -O0 level and -O1 level were different but from -O1 to -O3 gave same (wrong) results on gfortran, g77 and g95. I tested it with pgi fortran and got same (right) results. I checked gfortran 4.0.4, 4.1.2 and 4.2.0. I did not check gfortran 4.3. 2007/6/18, Tobias Burnus <burnus@net-b.de>: > Sunjoong Lee wrote: > > I had compiled a legacy fortran77 code and foud a bug; > > $ gfortran -o TMalign TMalign.f > > $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali > > Aligned length= 89, RMSD= 6.41, TM-score= 0.24257, ID=0.042 > > $ gfortran -O0 -o TMalign TMalign.f > > $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali > > Aligned length= 91, RMSD= 6.35, TM-score=0.24762 , ID=0.024 > I find this difference very odd as "-O0" is the default optimization for > gfortran. Other than that I get always the same result ("91") with all > -O levels I tried with gfortran 4.3, 4.2.0, 4.1.3 and ifort. > > Which version of the compiler are you using on which platform. (Use > "gfortran -v" to shows this information.) > Can you also show what "alias gfortran" (or "type gfortran") shows? Just > to make sure there is no alias which adds options. > > Tobias >
I did not make as many tests as Tobias, but it "works" for me on PPC with g77, xlf, g95, and gfortran.
Thanks again, Tobias $ uname -a Linux newton 2.6.12-gentoo-r10 #1 Sun Sep 4 13:29:37 KST 2005 i686 Intel(R) Pentium(R) 4 CPU 2.40GHz GenuineIntel GNU/Linux $ gfortran -v Using built-in specs. Target: i686-pc-linux-gnu Configured with: /var/tmp/portage/sys-devel/gcc-4.1.2/work/gcc-4.1.2/configure --prefix=/usr --bindir=/usr/i686-pc-linux-gnu/gcc-bin/4.1.2 --includedir=/usr/lib/gcc/i686-pc-linux-gnu/4.1.2/include --datadir=/usr/share/gcc-data/i686-pc-linux-gnu/4.1.2 --mandir=/usr/share/gcc-data/i686-pc-linux-gnu/4.1.2/man --infodir=/usr/share/gcc-data/i686-pc-linux-gnu/4.1.2/info --with-gxx-include-dir=/usr/lib/gcc/i686-pc-linux-gnu/4.1.2/include/g++-v4 --host=i686-pc-linux-gnu --build=i686-pc-linux-gnu --disable-altivec --enable-nls --without-included-gettext --with-system-zlib --disable-checking --disable-werror --enable-secureplt --disable-libunwind-exceptions --disable-multilib --enable-libmudflap --disable-libssp --disable-libgcj --enable-languages=c,c++,fortran --enable-shared --enable-threads=posix --enable-__cxa_atexit --enable-clocale=gnu Thread model: posix gcc version 4.1.2 (Gentoo 4.1.2) I had tested it with this step; $ gfortran -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ gfortran -O1 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ wget -O - http://ftp.g95.org/g95-x86-linux.tgz | tar xvfz - $ export PATH="$PATH:./g95-install/bin" $ g95 -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ g95 -O1 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ emerge =sys-devel/gcc-4.2.0 # install gcc-4.2.0 and gfortran $ gcc-config i686-pc-linux-gnu-4.2.0 $ source /etc/profile $ epm -e =sys-devel/gcc-4.1.2 # remove gcc-4.1.2 $ gfortran -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ gfortran -O1 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ emerge =sys-devel/gcc-4.0.4 # install gcc-4.0.4 and gfortran $ gcc-config i686-pc-linux-gnu-4.0.4 $ source /etc/profile $ epm -e =sys-devel/gcc-4.2.0 # remove gcc-4.2.0 $ gfortran -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ gfortran -O1 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali After then, I loged in another machine; $ uname -a Linux gene.kias.re.kr 2.6.3-gene #3 SMP Wed Jan 19 00:10:01 KST 2005 i686 unknown unknown GNU/Linux $ g77 -v Reading specs from /usr/lib/gcc-lib/i586-mandrake-linux-gnu/3.3.2/specs Configured with: ../configure --prefix=/usr --libdir=/usr/lib --with-slibdir=/lib --mandir=/usr/share/man --infodir=/usr/share/info --enable-shared --enable-threads=posix --disable-checking --enable-long-long --enable-__cxa_atexit --enable-clocale=gnu --enable-languages=c,c++,ada,f77,objc,java,pascal --host=i586-mandrake-linux-gnu --with-system-zlib Thread model: posix gcc version 3.3.2 (Mandrake Linux 10.0 3.3.2-6mdk) $ g77 -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ g77 -O1 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ pgf77 -O0 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali $ pgf77 -O1 -o TMalign TMalign.f $ ./TMalign 1aquA.pdb 1avaC.pdb | grep ^Ali Which do you think the reson for this difference be for fortran compiler or c compiler? 2007/6/19, Tobias Burnus <burnus@net-b.de>: > As said: It works here with 4.1.3 20070521, 4.2.1 20070604 and 20070618, > 4.3.0 20070618. It also work with my g95, Intel Fortran and sunf95 > compilers. In all cases I get: > > Aligned length= 91, RMSD= 6.35, TM-score=0.24762, ID=0.024 > > I'm on x86_64-unknown-linux-gnu (openSUSE Factory) and it works both in > 32 and 64 bit mode with no option, -O0, -O1, -O2, -O3, -O3 -ffast-math. > > > What is your platform (CPU type and operating system) and what is your > exact version of gfortran? (Both information is shown by "gfortran -v".) > (You are really only using "-O1" and not any other flags, are you?) > > Tobias >
> $ gfortran -v > Target: i686-pc-linux-gnu > gcc version 4.1.2 (Gentoo 4.1.2) I can reproduce this on my Pentium M (openSUSE 10.2) with Target: i586-suse-linux gcc version 4.1.2 20061115 (prerelease) (SUSE Linux) and with FX's Target: i386-pc-linux-gnu gcc version 4.3.0 20070618 (experimental) Interestingly, it works on x86_64-unknown-linux-gnu even with -m32. (I somehow fail to run FX's i386 compiler on x86-64.) I think this could be a middle-end or target problem as it is that target dependent. On the other hand, gfortran (and g95) could simply produce a wrong generic code.
I checked which part of TMalign.f make optimizaton wrong; In DP subroutine, DO j=1,NSEQ2 DO i=1,NSEQ1 D=VAL(i-1,j-1)+SCORE(i,j) H=VAL(i-1,j) if(DIR(i-1,j))H=H+GAP_OPEN V=VAL(i,j-1) if(DIR(i,j-1))V=V+GAP_OPEN IF((D.GE.H).AND.(D.GE.V)) THEN DIR(I,J)=.true. VAL(i,j)=D ELSE DIR(I,J)=.false. if(V.GE.H)then val(i,j)=v else val(i,j)=h end if ENDIF ENDDO ENDDO
I cut the bellow and made a new subroutine, then another part did not change on '-O0' and '-O1'; D=VAL(i-1,j-1)+SCORE(i,j) H=VAL(i-1,j) if(DIR(i-1,j))H=H+GAP_OPEN V=VAL(i,j-1) if(DIR(i,j-1))V=V+GAP_OPEN IF((D.GE.H).AND.(D.GE.V)) THEN DIR(I,J)=.true. VAL(i,j)=D ELSE DIR(I,J)=.false. if(V.GE.H)then val(i,j)=v else val(i,j)=h end if ENDIF
There are literal hundreds of warning given by ftnchek, and there appears to be an array bounds problem. troutmask:sgk[231] gfc4x -o z -O -fbounds-check TMalign.f troutmask:sgk[232] ./z 1aquA.pdb 1avaC.pdb | grep RMSD At line 2122 of file TMalign.f Fortran runtime error: Array reference out of bounds for array 'w', upper bound +of dimension 1 exceeded
Yes, I agree that program is not beautiful and I know the the array 'w' of 'u3b' subroutine problem; I think the author of u3b use w(1) as pointer. However, the wrong generation of optimized code occurs in 'DP' subroutine. The array of DP have fixed boundary. (In reply to comment #10) > There are literal hundreds of warning given by ftnchek, and there > appears to be an array bounds problem. > > troutmask:sgk[231] gfc4x -o z -O -fbounds-check TMalign.f > troutmask:sgk[232] ./z 1aquA.pdb 1avaC.pdb | grep RMSD > At line 2122 of file TMalign.f > Fortran runtime error: Array reference out of bounds for array 'w', upper bound > +of dimension 1 exceeded >
(In reply to comment #11) > Yes, I agree that program is not beautiful > and I know the the array 'w' of 'u3b' subroutine problem; > I think the author of u3b use w(1) as pointer. Change the 1 to *. > However, > the wrong generation of optimized code occurs in 'DP' subroutine. > The array of DP have fixed boundary. Do you get wrong results if you add the -ffloat-store option? Can you also try the -fdefault-real-8 option?
The '-ffloat-store' option works! Thank you. However that gave me some quenstions; Is that feature or bug? There is many floating point operations of course. Why the only one specific resion make problem? That is not set when optimize. So, what's the disadvantage and why the 'fortran' with many floating point operation dosn't take it when optimize? Why x86 need that option but x86_64 or PPC do not need? (In reply to comment #12) > Do you get wrong results if you add the -ffloat-store option? > Can you also try the -fdefault-real-8 option?
(In reply to comment #13) > The '-ffloat-store' option works! Thank you. > > However that gave me some quenstions; > > Is that feature or bug? It is a 'feature' of the i386 class of cpu. See PR 323 for details. > There is many floating point operations of course. > Why the only one specific resion make problem? The floating point operations are done with 80-bit precision and the intermediate results that are left in a register have the extra precision. Double precision has 53 bits of precision and single precision has 24 bits. The mixed mode arithmetic that ftnchek warns about can cause the type of problem you are seeing. The -ffloat-store option forces the contents in a register after a floating pointing operation to be written to main memory and then re-read. This removes the extra precision. > That is not set when optimize. > So, what's the disadvantage and > why the 'fortran' with many floating point operation > dosn't take it when optimize? You'll need to review the code. I'd suggest first eliminating the warns produced by ftnchek. If the problems disappear, then be happy. If the problems are still present, you'll need to review the floating point operations in the code. > Why x86 need that option but x86_64 or PPC do not need? These cpus do not store intermediate results in 80-bit register.
So this is just a dup of bug 323 so closing as such. *** This bug has been marked as a duplicate of 323 ***
Thank all of you. I could understand what make it different. There is no 'volatile' statement in fortran77 syntax of gfortran. Of course, volatile is not fortran77 standard, I think, but a certian implimentation support volatile. http://web.utk.edu/~prdaves/Computerhelp/Fortran_Reference/fortran_statements.htm I made a bellow c function and checked it happends. Yes, the problem is same but in the c function, I can use 'volatile' keyword and be happy. (I hope the next version of gfortran support volatile statement in fortran 77.) C language version of decide subroutine (or decide_ function); #define NMAX 5000 extern struct { float SCORE[NMAX][NMAX]; float GAP_OPEN; int INVMAP[NMAX]; } dpc_; void decide_(int *i, int *j, int iDIR[NMAX + 1][NMAX + 1], float VAL[NMAX + 1][NMAX + 1]) { volatile float D; float H,V; D = VAL[*j - 1][*i - 1] + dpc_.SCORE[*j - 1][*i - 1]; if(iDIR[*j][*i - 1] == 1) H = VAL[*j][*i - 1] + dpc_.GAP_OPEN; else H = VAL[*j][*i - 1]; if(iDIR[*j - 1][*i] == 1) V = VAL[*j - 1][*i] + dpc_.GAP_OPEN; else V = VAL[*j - 1][*i]; if((D >= H) && (D >= V)) { iDIR[*j][*i] = 1; VAL[*j][*i] = D; } else { iDIR[*j][*i] = 0; if(V >= H) VAL[*j][*i] = V; else VAL[*j][*i] = H; } } DP subroutine use above decide subroutine; SUBROUTINE DP(NSEQ1,NSEQ2) PARAMETER(nmax=5000) common/dpc/score(nmax,nmax),gap_open,invmap(nmax) dimension iDIR(0:nmax,0:nmax),VAL(0:nmax,0:nmax) REAL H,V C** initialize the matrix: val(0,0)=0 do i=1,nseq1 idir(i,0)=0 val(i,0)=0 enddo do j=1,nseq2 idir(0,j)=0 val(0,j)=0 invmap(j)=-1 enddo C** decide matrix and path: DO j=1,NSEQ2 DO i=1,NSEQ1 call decide(i,j,iDIR,VAL) ENDDO C** extract the alignment: i=NSEQ1 j=NSEQ2 DO WHILE((i.GT.0).AND.(j.GT.0)) IF(iDIR(i,j).eq.1)THEN invmap(j)=i i=i-1 j=j-1 ELSE H=VAL(i-1,j) if(iDIR(i-1,j).eq.1)H=H+GAP_OPEN V=VAL(i,j-1) if(iDIR(i,j-1).eq.1)V=V+GAP_OPEN IF(V.GE.H) THEN j=j-1 ELSE i=i-1 ENDIF ENDIF ENDDO return END
(In reply to comment #16) > Thank all of you. > I could understand what make it different. > > There is no 'volatile' statement in fortran77 syntax of gfortran. > Of course, volatile is not fortran77 standard, I think, > but a certian implimentation support volatile. You need to update the Fortran Standard that you use. Fortran 77 hasn't been the standard since about 1990. In fact, you're 3 standard been! There was Fortran 90, which was replaced by Fortran 95, which was replaced by Fortran 2003. Fortran 2003 has the VOLATILE attribute and VOLATILE statement. Guess what?? No, go ahead and guess! gfortran supports this feature. You need to go read Goldberg's paper about floating point arithmetic.
I appreciate kargl's comments; they were helpful. I had known there is VOLATILE attribute in new Fortran standard but I had worked with "LEGACY" fortran77 program! I'll write C code if I shuld write one; that is more compatable for me. However there are many legacy fortran 77 code in my field and I don't want to rewrite it. Small adjustment like '-ffloat-store' when compiling by gfortran instead of pgf77 is acceptable for I understood the behavior of register kargl told me. (I had read Goldberg's briefly but that is not point.) I had ONLY HOPEd VOLATILE statement in fortran 77 EXTENSION of gfortran. I thought that would be convenient on small modification of legacy fortran 77 program. (In reply to comment #17) > You need to update the Fortran Standard that you use. Fortran > 77 hasn't been the standard since about 1990. In fact, you're > 3 standard been! There was Fortran 90, which was replaced by > Fortran 95, which was replaced by Fortran 2003. Fortran 2003 > has the VOLATILE attribute and VOLATILE statement. Guess what?? > No, go ahead and guess! gfortran supports this feature. > > You need to go read Goldberg's paper about floating point arithmetic. >
(In reply to comment #18) > > I had ONLY HOPEd VOLATILE statement in fortran 77 EXTENSION of gfortran. > I thought that would be convenient > on small modification of legacy fortran 77 program. You've completed missed the point! gfortran is a Fortran 95 compiler with some Fortran 2003 extensions. A valid Fortran 77 code is a valid Fortran 95 code. You can add VOLATILE statements to your old legacy code and gfortran will compile it.
$ gfortran -O1 -o TMalign TMalign.f In file TMalign.f:2005 VOLATILE D 1 Error: Unclassifiable statement at (1) (In reply to comment #19) > You can add VOLATILE statements to your old legacy code and > gfortran will compile it.
(In reply to comment #20) > $ gfortran -O1 -o TMalign TMalign.f > In file TMalign.f:2005 > > VOLATILE D > 1 > Error: Unclassifiable statement at (1) Sigh. Try it with gfortran 4.3.0.
I checked VOLATILE statement passing in gfortran 4.3 . Unfortunately I had failed to buid gfortran (and gcc) 4.3 in my i686, I only checked that passing with ia64 binary on another ia64 machine. (There is no gfortran 4.3 binary on i686.) However, I think and hope it would work after my success of building on i686. (In reply to comment #21) > (In reply to comment #20) > > $ gfortran -O1 -o TMalign TMalign.f > > In file TMalign.f:2005 > > > > VOLATILE D > > 1 > > Error: Unclassifiable statement at (1) > > Sigh. Try it with gfortran 4.3.0.