[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: Possible loss of accuracy
From: |
Sergei Steshenko |
Subject: |
Re: Possible loss of accuracy |
Date: |
Wed, 15 May 2013 09:26:31 -0700 (PDT) |
--- On Wed, 5/15/13, Marco Caliari <address@hidden> wrote:
> From: Marco Caliari <address@hidden>
> Subject: Possible loss of accuracy
> To: address@hidden
> Date: Wednesday, May 15, 2013, 2:27 AM
> Dear all,
>
> if I compute
>
> t=50
> (1-2^(-t))^(2^t)
>
> I get (different versions, machines, ...)
>
> ans = 3.67879441128629e-01
>
> which is correct to the 10th digit (compared with
> www.wolframalpha.com). I can stay with it, but Matlab is
> correct to the 15th digit. While trying to understand, I
> wrote the following Fortran code
>
> program test
> integer t
> t = 50
> write(6,'(1x,e24.16)') (1.d0-2.d0**(-t))**(2.d0**t)
> write(6,'(1x,e24.16)') (1.d0-2.d0**(-50))**(2.d0**50)
> stop
> end
>
> If I compile with gfortran test.f (gfortran 4.3.3) and run,
> I get
>
> 0.3678794411286288E+00 (~ Octave value)
> 0.3678794411714422E+00 (~ Matlab value)
>
> and if I compile with gfortran -O2 test.f and run, I get
>
> 0.3678794411714422E+00
> 0.3678794411714422E+00
>
> So, I understand that compiler optimizations are involved.
> On the other hand, I compiled Octave with -O2 (and other
> flags).
>
> Can anyone point me to the right direction?
>
> Thanks,
>
> Marco
> _______________________________________________
> Help-octave mailing list
> address@hidden
> https://mailman.cae.wisc.edu/listinfo/help-octave
>
On my self-built octave-3.6.4 built with self-built gcc-4.5.4:
"
octave:5> output_precision(20)
octave:6> t=50; (1-2^(-t))^(2^t)
ans = 3.6787944112862880353e-01
octave:7>
".
And with self-built perl-5.16.3 built with the same self-built gcc-4.5.4:
"
address@hidden:~/junk> ~/AFSWD/install/perl-5.16.3/binsh/perl -e 'my $t=50;
print sprintf("%1.20e\n", ((1-2 ** (-$t)) ** (2 ** $t)))'
3.67879441128628803526e-01
".
So, Octave is in agreement with Perl, hence I think it's gcc (and possibly
glibc) issues.
The machine is a 32 bit one.
Regards,
Sergei.