forums.silverfrost.com Forum Index forums.silverfrost.com
Welcome to the Silverfrost forums
 
 FAQFAQ   SearchSearch   MemberlistMemberlist   UsergroupsUsergroups   RegisterRegister 
 ProfileProfile   Log in to check your private messagesLog in to check your private messages   Log inLog in 

Different results for 32-Bit and 64-Bit versions

 
Post new topic   Reply to topic    forums.silverfrost.com Forum Index -> Support
View previous topic :: View next topic  
Author Message
Theo1002



Joined: 19 Apr 2021
Posts: 15
Location: Ingelheim am Rhein, Germany

PostPosted: Wed Apr 21, 2021 10:27 pm    Post subject: Different results for 32-Bit and 64-Bit versions Reply with quote

Using Silverfrost ftn95 for reactivating older Fortran programs which run in the 80s on an IBM 3032 Mainframe generally works fine, when replacing IBM specific commands/conventions.
Currently, however, I'm stuck with a program calculating matrix elements. What puzzles me is the fact, that the 32Bit-Version delivers different results compared with the results of the 64Bit-Version (I'm using Plato for compilation and SDBG version 8.70.00).
Whereas the 32Bit-version produces correct results, the ones calculated with the 64Bit version are obvioulsy wrong. Trying to isolate the problem, kept me busy the last days, but is not so easy due to 20 subroutines and many function calls. Since I didn't succeed to isolate the root cause of the differences between the results of the 32- and 64-Bit versions up to now, I would like to ask the following questions:
Did anybody experience similar deviations between results of 32- and 64-Bit versions?
Can certain compiler options help to isolate the problem?
Please note, that I have no experience with command line use of SDBG, but use the Plato environment.
@mecej4: Thank you for the response and information contained therein concerning possible reasons for the different results. As suggested I add the part of source code here, which yields the different results.
Update: Meanwhile I could isolate the problem to origin from the subroutine wavefun. When calling the Whittaker function, variable GAMBA differs (as result from calling GAMCAR (Gamma function). The input used for WAVEFUN is as follows:
    -1.500D+00 -1 1.840D+02 3.286D-01 3.002D-01 5.000D+03 0 0.000D+0
    3.2858593764968097D-01 3.0015179886716270D-01

The subroutine WAVEFUN (with many writes to file 14 for test purposes) is listed below:

[code:1:aabccb7571]C ******************************************************************
C ******************************************************************
C SUBROUTINE WAVFUN(E,KAPPA,Z,U1,U2,R0,ITEST,DELTA)
C ******************************************************************
C ******************************************************************
C THIS PROGRAM COMPUTES THE RADIAL KONTINUUM WAVEFUNCTIONS OF THE
C DIRAC EQUATION FOR EACH Z WITH RESPECT OF AN EXTENDED NUCLEUS.
C U1 AND U2 ARE THE RADIAL KONTINUUMWAVEFUNCTIONS.
C SEE FOR EXAMPLE EQ. 5.4' OF REF 2.
C THE RADIUS OF THE NUCLEUS IS GIVEN BY R=R0*A**(1./3.)
C UNITS ARE H/(2.*PI) = C = M = 1.
C RF IS R IN FM.
C REF 1. B. MUELLER, J. RAFELSKI AND W. GREINER
C INSTITUT FUER THEORETISCHE PHYSIK - FRANKFURT
C PREPRINT (1973)
C REF 2. E. M. ROSE, RELATIVISTISCHE ELEKTRONENTHEORIE.
C REF 3. B. MUELLER, J. RAFELSKI, AND W. GREINER
C Z. PHYSIK 257,183 (1972)
IMPLICIT REAL*8 (A-H,O-Z)
COMPLEX*16 PHI1,PHI2,FAK3,FAK4,GAM,CDSQRT,ABC,FAK1,FAK2,DCONJG
COMMON /ZELLE/A1,A2
OPEN (5,file='D:\FORTRAN\MEDDAT\TST\FT05WF.DAT',status='READONLY')
OPEN (6,file='D:\FORTRAN\MEDDAT\TST\FT06WF.DAT',status='OLD')
OPEN (14,file='D:\FORTRAN\MEDDAT\TST\FT14WF.DAT',status='OLD')
READ (5,*) E,KAPPA,Z,U1,U2,R0,ITEST,DELTA,A1,A2
WRITE (14,1680) E,KAPPA,Z,U1,U2,R0,ITEST,DELTA
1680 FORMAT (1X,'Entering WAVEFUN, E,KAPPA,Z,U1,U2,R0,ITEST,DELTA',
&1X,1PD10.3,I3/1X,1P4D15.3,I6,1PD10.3)
COMWAV=386.1592D0
ALPHA=1./137.03602D0
PI=3.1415926535897D0
Back to top
View user's profile Send private message
Theo1002



Joined: 19 Apr 2021
Posts: 15
Location: Ingelheim am Rhein, Germany

PostPosted: Wed Apr 21, 2021 10:41 pm    Post subject: Source code Reply with quote

The source code is located unter the following link:
https://www.magentacloud.de/share/gnh315jfie
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1580

PostPosted: Wed Apr 21, 2021 11:53 pm    Post subject: Reply with quote

Theo, you posted only the source file and one data file to MagentaCloud. The file OPEN statements are inconsistent with the usage of the files. I corrected this by changing 'OLD' to 'REPLACE' in the OPEN statements for units 6 and 14.

As you may have noticed, there is a limit of about 100 lines per post, so long program sources may lose major portions when posted inline. For short code extracts, please select the source lines with the mouse and click on the Code button. The lines will then be formatted nicely.

The function GAMCAR takes one argument, CMX. This variable is both read and updated in the function. This is usually a cause of trouble, because the normal expected behaviour of a function is to return a function value without changing any of its arguments. In your program, this function GAMCAR is also referenced with an expression as an argument (on line 419). With most modern compilers, you are not allowed to place a value into an expression -- you can only place values into variables. I corrected this by changing the name of the dummy argument of GAMCAR to CMXARG, and setting CMX=CMXARG at the beginning of the function.

After these changes, I compiled and ran the program. The program printed "STOP: 111" and stopped, since the input data provides '0' for ITEST.

The results in the output file generated, FT14WF.DAT, was nearly the same with two different compilers: FTN95 8.71, 64-bit, and Lahey LF71, 32-bit. There were small differences in the least significant digits, but that is not unusual with floating point calculations.

I then compared the output of 32-bit FTN95 and 64-bit FTN95. Again, the results agreed to at least 12 significant digits.

What next? What did you do to produce the different results that you did not feel were correct?
Back to top
View user's profile Send private message
Theo1002



Joined: 19 Apr 2021
Posts: 15
Location: Ingelheim am Rhein, Germany

PostPosted: Thu Apr 22, 2021 10:50 am    Post subject: Reply with quote

mecej4, first of all thank you for your second timely response and reviewing the source code! I added the two files delivering different results when compiling the original source code with 32 Bit, resp. 64 Bit compiler (FT14WF-32BIT.DAT and FT14WF-64BIT.DAT). The differences in GAMBA are as follows (first printout):

(1.47646338812876D-02,-1.63783687568201D-02) 32-Bit
(7.38231694064379D-03,-8.18918437841007D-03) 64-Bit

This answers your last question. By error propagation from this subroutine I got faulty results in the matrix elements.

I was'nt aware that passing arguments forth and back via the subroutine parameter list is a problem for current compilers, since that was common use (and allowed) in IBM- Fortran (about 40 years ago..). So thanks a lot for this valuable hint! Same with line 419: only argument, no expression in function call.

After changing this according to your recommendation, I got agreeing results - within 11-12 significant digits (c.f. file FT14WF-64Bit-new.DAT) between FTN95 32- and 64- Bit runs. Stop 111 is a regular RETURN in the subroutine.
So the FTN95 32Bit compiler obviously tolerates the mentioned coding deficiencies, wheras the 64Bit compiler doesn't.

I guess I wouldn't have figured it out without your help !
What next? I will search the complete source code concerning modification of passed arguments in subroutine calls (hours of enjoyment;-) plus possible expressions in function argument lists, correct and than check the final result.
Thank you very much for your help, mecej4!

Kind regards
Theo
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1580

PostPosted: Thu Apr 22, 2021 11:07 am    Post subject: Reply with quote

Theo, I am gratified to see you catching up so quickly, and I have good news for you. Silverfrost FTN95 is very good at catching bugs in Fortran programs at compile time (i.e., before even running the program), but there are many bugs that can only be caught at run time, and FTN95 provides facilities to catch those as well.

Try the /debug, /check and /undef options. When your code contains more Fortran 95 /2003 features, you may also find /checkmate useful.

Here is a buggy but short program that may interest you. Try running it with different options, 32 and 64 bit, and you may find some of the outcomes quite puzzling and surprising. In this connection, see also https://medium.com/@patkerr/should-the-value-of-pi-change-f18f41599919 .

Code:
program OneIsTwo
implicit none
integer One, Two
One = 1
Two = 2
call Sub(One, Two)
Print *, 'One = ',One, 'Two = ',Two
call Sub(1, 2)
Print *, 'Later, One = ',1, ' and Two = ',2
end program

subroutine sub(A,B)
implicit none
integer A,B,C
C=A
A=B
B=C
return
end subroutine


What would the old IBM computer output for this program (after reformatting the source to Fortran 4 or Fortran 77 conventions to make it compatible)?
Back to top
View user's profile Send private message
wahorger



Joined: 13 Oct 2014
Posts: 1011
Location: Morrison, CO, USA

PostPosted: Thu Apr 22, 2021 5:48 pm    Post subject: Reply with quote

Theo, is it just serendipity that the answers differ by a factor of 2?

Bill
Back to top
View user's profile Send private message Visit poster's website
Theo1002



Joined: 19 Apr 2021
Posts: 15
Location: Ingelheim am Rhein, Germany

PostPosted: Thu Apr 22, 2021 6:00 pm    Post subject: Reply with quote

Well the first result is as expected in good old IBM Fortran (One = 2 and Two =1) and the second yields an error - as expected.

Using checkmate indeed delivers more information:
    Runtime error from program:d:\fortran\me-dive\fixedformat1.exe
    Run-time Error
    *** Error 14, Attempt to alter an actual argument that is a constant, an expression, an alias, an INTENT(IN) argument, or a DO variable

    SUB - in file fixedformat1.for at line 16 [+009f]

    ONEISTWO - in file fixedformat1.for at line 8 [+01af]

I have used the /debug, /check and /undef parameters and eliminated tons of undefined variables in the FT77 pgms (IBM Fortran used to initialize every variable nicely with 0), but that didn't give me a hint to the problem with the argument list in the subroutine calls.

Thanks again for your help - I have truncated the source code in magenta cloud to the relevant functions.

Kind regards
Theo
Back to top
View user's profile Send private message
Theo1002



Joined: 19 Apr 2021
Posts: 15
Location: Ingelheim am Rhein, Germany

PostPosted: Fri Apr 23, 2021 2:38 pm    Post subject: Reply with quote

Good point Bill (eagle-eye)! I was so busy chasing the occurence, where the differences between the two compiler versions occur, that I didn't even notice that the results differ by a factor 2.

The factor 2 indeed seems to be serendipity, since (to my knowledge) there is no recursion relation for the Gamma functions in the form Gamma(B-A)=Gamma(B)/N or Gamma(A)/N [in case GAMCAR(B-A) would be interpreted either as GAMCAR(B) or GAMCAR(A)] for complex arguments A,B and integer N. The difference origins from the faulty argument passing at FUNCTION GAMCAR(CMX) via CMX as described by mecej4. It seems that this leads to returning a faulty function value by a factor 2.

Best regards
Theo
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1580

PostPosted: Fri Apr 23, 2021 4:51 pm    Post subject: Reply with quote

Theo and Bill,

That factor of 2 that you noticed is just happenstance. When a Fortran program is in violation of the standard, most often the compiler and RTL may do anything, including

    1. Refuse to compile the program.
    2. Produce an EXE that runs a bit, gives an error message and stops.
    3. Produce an EXE that hangs up, and you may have to use the OS tools to abort it.
    4. Produce incorrect results. If they are only slightly incorrect, the user may not even notice that something is amiss.
    5. Produce correct results.
    6. Start WW-III.

and so on.

The Fortran standard simply uses the phrase "behavior is undefined" to cover all these possibilities.


Last edited by mecej4 on Fri Apr 23, 2021 5:16 pm; edited 1 time in total
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1580

PostPosted: Fri Apr 23, 2021 4:58 pm    Post subject: Function GAMCAR in modern Fortran Reply with quote

Theo's IBM 3033 code contains a surprisingly good and efficient function subprogram for computing the Gamma function of a complex argument within certain ranges of the argument. I found it a striking instance where the facilities of modern Fortran can enable the programmer to express the algorithm in a compact, efficient and readable form, in contrast to the "spaghetti" style of the original.

Code:
!     +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
      complex(dble_complex) function gamcar (cmx)
      USE kind_param, ONLY:  dble_complex, double
!     +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
!     THIS SUBPROGRAM EVALUATES THE GAMMA FUNCTION FOR ANY
!     COMPLEX ARGUMENT BY USE OF AN ASYMPTOTIC EXPANSION.
!     THIS PROGRAM HAS BEEN CHECKED FOR 1. LE RE(Z) LE 2. AND
!     0.1 LE IM(Z) LE 10.0.
!     ERROR IS LESS THAN 10**-6
!**************************************************************
      implicit none
      integer :: j, k
      real(double) :: t
      complex(dble_complex) :: cmx, comx, tx, gam1
!-----------------------------------------------
      k = nint(dble(cmx))
      comx = cmx + k
      tx   = comx + 5.5
      gam1 = exp((comx + 0.5)*log(tx) - tx)*2.50662827465D0
      gamcar = (1. + 76.18009173D0/(comx + 1.) - 86.50532033D0/(comx + 2.) &
                   + 24.01409822D0/(comx + 3.) - 1.231739516D0/(comx + 4.) &
                   + 0.12085003D-2/(comx + 5.) -   0.536382D-5/(comx + 6.) &
               )*gam1/product([(cmx+j, j=0,k)])
      return
      end function gamcar


Contrast this with the original (with potentially illegal modification of function argument), continued in the next post in the thread (forum limit)

[CONTINUED]


Last edited by mecej4 on Fri Apr 23, 2021 8:32 pm; edited 1 time in total
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1580

PostPosted: Fri Apr 23, 2021 5:00 pm    Post subject: Reply with quote

[CONTINUED]
Code:
C     +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
      FUNCTION GAMCAR(CMX)                                             
C     +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
C     THIS SUBPROGRAM EVALUATES THE GAMMA FUNCTION FOR ANY             
C     COMPLEX ARGUMENT BY USE OF AN ASYMPTOTIC EXPANSION.               
C     THIS PROGRAM HAS BEEN CHECKED FOR 1. LE RE(Z) LE 2. AND           
C     0.1 LE IM(Z) LE 10.0.                                             
C     ERROR IS LESS THAN 10**-6                                         
      IMPLICIT REAL*8 (A-H,O-Z)                                         
      COMPLEX*16 GAMCAR                                                 
      COMPLEX*16 CMX,COMX,GAMX,E1,E11,TX,GAM1,SUMME,FAK       
      SUMME=CMX                                                           
      FAK=(1.,0.)                                                       
      IF (DBLE(CMX).LE.0.) GO TO 5                                     
   12 CONTINUE                                                         
      COMX=CMX-1.                                                       
      E1=COMX+.5                                                       
      TX=COMX+5.5                                                       
      E11=E1*LOG(TX)                                                 
      GAM1=EXP(E11-TX)*2.50662827465D0                               
      GAMX=GAM1*(1.+76.18009173D0/(COMX+1.)-86.50532033D0/(COMX+2.)     
     &+24.01409822D0/(COMX+3.)-1.231739516D0/(COMX+4.)+.12085003D-2/(COM
     &X+5.)-.536382D-5/(COMX+6.))                                       
      GAMCAR=GAMX                                                       
      GAMCAR=GAMCAR*FAK                                                 
      GO TO 2                                                           
    5 CONTINUE                                                         
      J=0                                                               
    6 CMX=CMX+1.                                                       
      J=J+1                                                             
      T=DBLE(CMX)                                                     
      IF (T.LE.0) GO TO 6                                               
      DO 11 K=1,J                                                       
      FAK=FAK*1./(SUMME+K-1.)                                     
   11 CONTINUE                                                         
      GO TO 12                                                         
    2 CMX=SUMME
      RETURN                                                           
      END
Back to top
View user's profile Send private message
Theo1002



Joined: 19 Apr 2021
Posts: 15
Location: Ingelheim am Rhein, Germany

PostPosted: Sat Apr 24, 2021 10:32 am    Post subject: Reply with quote

mecej4, thank you for explaining what happens when a Fortran program is in violation of standards.

The specific code for evaluating the Gamma function is a module which was not coded by myself but by one of the members of the Theoretical Physics group in Frankfurt (larger programs were jointly coded by the members of the Atomic Physics group). At that time I would have code it simlarly within FT77 Standard.

Your recoding is fascinating short and elegant ! (maybe I should have a deeper look into FTN95 coding standards). I will definitly test and use this.

Best regards
Theo
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1580

PostPosted: Sun Apr 25, 2021 1:14 pm    Post subject: Reply with quote

The IMSL library has the generic function GAMMA -- not surprising --, so I compared your results with what IMSL gives. If you have IMSL or have access to a computer with IMSL, you can test this code:

Code:
PROGRAM TGAM
USE GAMMA_INT
IMPLICIT NONE
integer, parameter :: double = kind(0d0)
!
complex(double) :: ab,gab
complex sab,gimsl
ab = (0d0,2.69747d0)
sab = ab
gab = gamcar(ab)
gimsl = gamma(sab)
print 10,ab,gab,gimsl
10 format(/,' Argument : ',2ES13.5,/,' GamCar   : ',2ES13.5,/, &
            ' IMSL     : ',2ES13.5)

CONTAINS
!     +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
      complex(double) function gamcar (cmx)
!     +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-
!     THIS SUBPROGRAM EVALUATES THE GAMMA FUNCTION FOR ANY
!     COMPLEX ARGUMENT BY USE OF AN ASYMPTOTIC EXPANSION.
!     THIS PROGRAM HAS BEEN CHECKED FOR 1. LE RE(Z) LE 2. AND
!     0.1 LE IM(Z) LE 10.0.
!     ERROR IS LESS THAN 10**-6
!**************************************************************
      implicit none
      integer :: j, k
      real(double) :: t, G0 = 2.50662827465D0
      real(double) :: G(6) = [76.18009173D0, -86.50532033D0, 24.01409822D0, &
                              -1.231739516D0, 0.12085003D-2, -0.536382D-5]
      complex(double) :: cmx, comx, tx, gam1
!-----------------------------------------------
      k = nint(dble(cmx))
      comx = cmx + k
      tx   = comx + 5.5
      gam1 = exp((comx + 0.5)*log(tx) - tx)*G0
      gamcar = (1. + sum([(G(j)/(comx + j), j = 1, 6)])) &
              *gam1/product([(cmx+j, j=0,k)])
      return
      end function gamcar
end program


Output:

Code:
T:\LANG\THEO\F90>f95 %f90flags% gamcar.f90 %link_fnl%

T:\LANG\THEO\F90>gamcar

 Argument :   0.00000E+00  2.69747E+00
 GamCar   :   1.47647E-02 -1.63785E-02
 IMSL     :   1.47647E-02 -1.63785E-02


For the argument that I used, at least, it appears that COMPLEX*16 may not be needed.
Back to top
View user's profile Send private message
Display posts from previous:   
Post new topic   Reply to topic    forums.silverfrost.com Forum Index -> Support All times are GMT + 1 Hour
Page 1 of 1

 
Jump to:  
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum


Powered by phpBB © 2001, 2005 phpBB Group