This is the mail archive of the cygwin mailing list for the Cygwin project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

Re: cygwin gcc: Different numerical results in thread vs in main()


Brian Dessent wrote:

Angelo Graziosi wrote:

... and in Fortran?

As long as you're using a recent gcc you can just use -mpc64.



How recent?


With GFortran 4.3.1 and

$ cat test_case.0.f90
program test_case
  implicit none
  integer :: k
  integer, parameter :: DP = kind(1.D0),&
       N = 29
  real(DP), parameter :: A(0:N) = &
       (/1.0000000000000D0,3.8860009363293D0,7.4167881843083D0,&
         9.8599837415463D0,10.431383465276D0,9.4077161304727D0,&
         7.5332956276775D0,5.4995844630326D0,3.7275104474917D0,&
         2.3766149078263D0,1.4394444941337D0,0.8344437543616D0,&
         0.4657058294147D0,0.2514185849207D0,0.1317920740387D0,&
         0.0672995854816D0,0.0335554996094D0,0.0163785206816D0,&
         0.0078368645385D0,0.0036846299276D0,0.0017020281304D0,&
         0.0007756352209D0,0.0003469750783D0,0.0001544224626D0,&
         0.0000666379875D0,0.0000295655909D0,0.0000118821415D0,&
         0.0000057747681D0,0.0000017502652D0,0.0000014682034D0/)
  real(DP), parameter :: BB(0:N-1) = &
       (/(((k+2)*A(k+2)-A(k+1)*(A(1)+(k+1))),k = 0,N-2),&
         (-A(N)*(A(1)+N))/)
  real(DP) :: b(0:N-1),bbb(0:N-1)
  do k = 0,N-2
     b(k) = (k+2)*A(k+2)-A(k+1)*(A(1)+(k+1))
  enddo
  b(N-1) = -A(N)*(A(1)+N)
  bbb(0:N-1) = (/(((k+2)*A(k+2)-A(k+1)*(A(1)+(k+1))),k = 0,N-2),&
       (-A(N)*(A(1)+N))/)
  do k = 0,N-1
     print *, k,b(k)-BB(k),b(k)-bbb(k),BB(k)-bbb(k)
  enddo
end program test_case

$ gfortran -mpc64 test_case.0.f90 -o test_case.0

$ ./test_case.0
0 -1.77635683940025046E-015 -1.77635683940025046E-015 0.0000000000000000
1 -3.55271367880050093E-015 -3.55271367880050093E-015 0.0000000000000000
2 -3.55271367880050093E-015 -3.55271367880050093E-015 0.0000000000000000
3 -7.10542735760100186E-015 -7.10542735760100186E-015 0.0000000000000000
4 0.0000000000000000 0.0000000000000000 0.0000000000000000
5 0.0000000000000000 0.0000000000000000 0.0000000000000000
...


Only 13 (over 28) 0. 0. 0.!


Angelo


--
Unsubscribe info:      http://cygwin.com/ml/#unsubscribe-simple
Problem reports:       http://cygwin.com/problems.html
Documentation:         http://cygwin.com/docs.html
FAQ:                   http://cygwin.com/faq/


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]