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 

Code from source with SPREAD causes memory corruption

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



Joined: 31 Oct 2006
Posts: 1064

PostPosted: Tue Dec 18, 2018 7:23 pm    Post subject: Code from source with SPREAD causes memory corruption Reply with quote

The test code below arose from an investigation of bugs in the venerable Eispack package at Netlib (test CHTEST). The old code has some undefined variables, and some subroutine calls with aliased arguments. FTN95 helped with the undefined variables, but not with the aliased arguments.

Refactored code from the investigation showed up an FTN95 bug. The bug is present in the 32 and 64 bit Version 8.30.290 compiler, but not in the 7.20 compiler. To keep the test code short, I have placed the data in an unformatted file (696 bytes), which can be downloaded from Dropbox: https://www.dropbox.com/s/i0loyu5ip2ya02q/sbug.sbn?dl=0 .

The code:
Code:
      program sbug
      implicit none
      integer, parameter :: double = kind(0.0d0)
      real(double) :: ar(5,5),ai(5,5),tau(2,5),zr(5,5),zi(5,5)
      integer :: n
      open(20,file='sbug.sbn',form='unformatted',status='old')
      read(20)n
      if(n /= 5)stop 'Incorrect n for data file'
      read(20)ar,ai,tau,zr
      close(20)
      call htribk (n, ar, ai, tau, zr, zi)
      end

      subroutine htribk(n, ar, ai, tau, zr, zi)
      implicit none
      integer, parameter :: double = kind(0.0d0)
      integer , intent(in) :: n
      real(double) , intent(in) :: ar(n,n), ai(n,n), tau(2,n)
      real(double) , intent(in out) :: zr(n,n)
      real(double) , intent(out)    :: zi(n,n)
      integer :: i, j, k
      real(double) :: h
      real(kind(0.0d0)) :: s1u(n), si1u(n)
!
      zi(:n,:) = -zr(:n,:)*spread(tau(2,:),dim = 2,ncopies = n)
      zr(:n,:) =  zr(:n,:)*spread(tau(1,:),dim = 2,ncopies = n)
      i = 2
      k = i - 1
      h = ai(i,i)
      s1u = sum(spread(ar(i,:k),dim = 2,ncopies = n)*zr(:k,:)- &
                spread(ai(i,:k),dim = 2,ncopies = n)*zi(:k,:), dim=1)
      si1u = sum(spread(ar(i,:k),dim = 2,ncopies = n)*zi(:k,:)+&
                spread(ai(i,:k),dim = 2,ncopies = n)*zr(:k,:), dim=1)
      s1u  = (s1u/h)/h
      si1u = (si1u/h)/h
      PRINT *,' i      S1u       Si1u:'
      print '(1x,i2,2ES13.4)',(j,s1u(j),si1u(j),j=1,n)
      stop
      end subroutine htribk


The correct output results, which FTN95 7.20 gives (and also Gfortran, etc., provided a suitable unformatted data file is used):
Code:
  i      S1u       Si1u:
  1   2.3405E-01   0.0000E+00
  2   8.5779E-02   0.0000E+00
  3   1.3610E-01   0.0000E+00
  4   1.0824E-01   0.0000E+00
  5   8.4260E-04   0.0000E+00

The results from FTN95 8.30.279:
Code:
32-bit

Default, /lgo
  i      S1u       Si1u:
  1   0.0000E+00   0.0000E+00
  2   9.8118-294   9.8089-294
  3   0.0000E+00   0.0000E+00
  4   9.8118-294   9.8089-294
  5   0.0000E+00   0.0000E+00

/checkmate /lgo
     Invalid floating point operation

/checkmate
  1   9.8040-294   9.8033-294
  2  -4.9760+234  -4.9760+234
  3   0.0000E+00   0.0000E+00
  4   9.8107-294   9.8078-294
  5   9.8040-294   9.8033-294

/opt
  i      S1u       Si1u:
  1   2.3405E-01   0.0000E+00
  2   8.5779E-02   0.0000E+00
  3   1.3610E-01   0.0000E+00
  4   1.0824E-01   0.0000E+00
  5   8.4260E-04   0.0000E+00      (all correct!)

[Message limit reached, to be continued in next post in thread]


Last edited by mecej4 on Tue Dec 18, 2018 7:59 pm; edited 1 time in total
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1064

PostPosted: Tue Dec 18, 2018 7:25 pm    Post subject: (continued from preceding post) Reply with quote

...CONTINUATION...

Code:
64-bit results

Default, /lgo
  i      S1u       Si1u:
  1   2.3405E-01   0.0000E+00      (correct)
  2   8.5779E-02   0.0000E+00         "
  3   1.3610E-01   0.0000E+00         "
  4   1.0824E-01   0.0000E+00         "
  5   2.3405E-01   0.0000E+00      (incorrect)

/checkmate
  same as for Default

/opt
  i      S1u       Si1u:
  1   2.3405E-01   0.0000E+00
  2   8.5779E-02   0.0000E+00
  3   1.3610E-01   0.0000E+00
  4   1.0824E-01   0.0000E+00
  5   8.4260E-04   0.0000E+00      (all correct)


NOTE: I rarely use the SPREAD intrinsic when writing code myself, and I noted that replacing the four statements containing SPREAD with the corresponding DO loops or other array expressions gives correctly functioning code. The SPREAD calls were introduced by an automatic Fortran clean-up utility; the original code was in Fortran 66.
Back to top
View user's profile Send private message
LitusSaxonicum



Joined: 23 Aug 2005
Posts: 1924
Location: Yateley, Hants, UK

PostPosted: Tue Dec 18, 2018 8:13 pm    Post subject: Reply with quote

It's a glitch in The Matrix, surely, when one gets that creepy sense of Déjà vu about a post in The Forum.
Back to top
View user's profile Send private message
mecej4



Joined: 31 Oct 2006
Posts: 1064

PostPosted: Tue Dec 18, 2018 8:21 pm    Post subject: Reply with quote

If you are interested in knowing how this got started, see the following thread by Jack Wolfman on Comp. Lang. Fortran:

https://groups.google.com/forum/#!topic/comp.lang.fortran/trzQK_bMipY

In the past, I had filed a couple of bug reports involving SPREAD. One of them involved the compiler rejecting valid code, and another a performance hit as a result of using SPREAD. I do not remember having seen wrong results caused by its use.

What is quite a surprise with the current bug report is that one gets correct results (mostly) with /OPT and really bad results with /CHECKMATE.
Back to top
View user's profile Send private message
PaulLaidler
Site Admin


Joined: 21 Feb 2005
Posts: 5713
Location: Salford, UK

PostPosted: Wed Dec 19, 2018 9:22 am    Post subject: Reply with quote

mecej4

Thank you for the feedback.

There are similarities to earlier outstanding bug reports and hopefully when they are fixed this one will also come out in the wash.
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