getting same output as gfortran, long_float

getting same output as gfortran, long_float

[Nasser M. Abbasi]

I am learning Ada again, and was wondering what do I need to change
in this small code to make it output the same precision as gfortran?

First I show the code, then explain more what the issue.

Fortran:

--------------------
PROGRAM foo
IMPLICIT NONE
INTEGER, PARAMETER :: DBL = KIND(1.D0)
REAL(KIND = DBL) :: x
   x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
PRINT *, x
END PROGRAM
-------------------

Now, gfortran has an option called -fdefault-real-8 which I thought
will map to Ada long_float.

https://gcc.gnu.org/onlinedocs/gfortran/Fortran-Dialect-Options.html
"-fdefault-real-8 Set the default real type to an 8 byte wide type."

When compiling and running the above Fortran program, this is the output

-------------------
gfortran -Wall -fdefault-real-8  foo2.f90
./a.out
    1.82898948331047096479195244627343369E-0003

gfortran -Wall  foo2.f90
./a.out
    1.8289894833104711E-003
----------------------------------------

It is the first result above I'd like to get with Ada. Here is
what I tried:

--------------------
with ada.text_io; use ada.text_io;
procedure foo_2 is
   x : Long_Float;
begin
   x := 12.0 * 0.0001/(1.0 * (1.0 - 0.1)**4 );
   put_line(long_float'image(x));
end foo_2;
----------------------

compile and run:

gnatmake foo_2.adb
gcc-4.8 -c foo_2.adb
gnatbind -x foo_2.ali
gnatlink foo_2.ali
./foo_2
  1.82898948331047E-03


side-by-side:
=============
Ada (long)                1.82898948331047E-03
Ada (long long)           1.82898948331047096E-03
Fortran:                  1.8289894833104711E-003
Fortran -fdefault-real-8: 1.82898948331047096479195244627343369E-0003

How to obtain the last result above in Ada?

ps. here is the long_long_float Ada version.

-------------
with ada.text_io; use ada.text_io;
procedure foo_2 is
   x : Long_Long_Float;
begin
   x := 12.0 * 0.0001/(1.0 * (1.0 - 0.1)**4 );
   put_line(Long_Long_Float'image(x));
end foo_2;
----------------

--Nasser

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Thu, 30 Apr 2015 16:17:49 -0500, Nasser M. Abbasi wrote:

> I am learning Ada again, and was wondering what do I need to change
> in this small code to make it output the same precision as gfortran?
> 
> First I show the code, then explain more what the issue.
> 
> Fortran:
> 
> --------------------
> PROGRAM foo
> IMPLICIT NONE
> INTEGER, PARAMETER :: DBL = KIND(1.D0)
> REAL(KIND = DBL) :: x
>    x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
> PRINT *, x
> END PROGRAM
> -------------------
> 
> Now, gfortran has an option called -fdefault-real-8 which I thought
> will map to Ada long_float.
> 
> https://gcc.gnu.org/onlinedocs/gfortran/Fortran-Dialect-Options.html
> "-fdefault-real-8 Set the default real type to an 8 byte wide type."
> 
> When compiling and running the above Fortran program, this is the output
> 
> -------------------
> gfortran -Wall -fdefault-real-8  foo2.f90
> ./a.out
>     1.82898948331047096479195244627343369E-0003
> 
> gfortran -Wall  foo2.f90
> ./a.out
>     1.8289894833104711E-003
> ----------------------------------------
> 
> It is the first result above I'd like to get with Ada. Here is
> what I tried:
> 
> --------------------
> with ada.text_io; use ada.text_io;
> procedure foo_2 is
>    x : Long_Float;
> begin
>    x := 12.0 * 0.0001/(1.0 * (1.0 - 0.1)**4 );
>    put_line(long_float'image(x));
> end foo_2;
> ----------------------
> 
> compile and run:
> 
> gnatmake foo_2.adb
> gcc-4.8 -c foo_2.adb
> gnatbind -x foo_2.ali
> gnatlink foo_2.ali
> ./foo_2
>   1.82898948331047E-03
> 
> 
> side-by-side:
> =============
> Ada (long)                1.82898948331047E-03
> Ada (long long)           1.82898948331047096E-03
> Fortran:                  1.8289894833104711E-003
> Fortran -fdefault-real-8: 1.82898948331047096479195244627343369E-0003
> 
> How to obtain the last result above in Ada?
>
> ps. here is the long_long_float Ada version.
> 
> -------------
> with ada.text_io; use ada.text_io;
> procedure foo_2 is
>    x : Long_Long_Float;
> begin
>    x := 12.0 * 0.0001/(1.0 * (1.0 - 0.1)**4 );
>    put_line(Long_Long_Float'image(x));
> end foo_2;
> ----------------

I don't know what FORTRAN does, but assuming, you are using Intel, the
longest floating-point available there is not, as you seem assume, of 16
bytes (quadruple-precision), but so-called extended precision, which is not
that long (10 bytes, namely). 

AFAIK, GNAT on x86 maps:

Float to single precision float (4 bytes)
Long_Float to double precision float (8 bytes)
Long_Long_Float to extended precision float (10 bytes, effectively)

You can try this program to see how many binary mantissa bits are available
for each type:
-----------------------------
with Text_IO; use Text_IO;
procedure Test is
begin
   Put_Line ("Float mantissa:" & Integer'Image (Float'Machine_Mantissa));
   Put_Line ("Long Float mantissa:" & Integer'Image
(Long_Float'Machine_Mantissa));
   Put_Line ("Long Long Float mantissa:" & Integer'Image
(Long_Long_Float'Machine_Mantissa));
end Test;
--------------------------
On an Intel machine it will print:

Float mantissa: 24
Long_Float mantissa: 53
Long_Long_Float mantissa: 64

If FORTRAN supports longer floating point numbers (e.g.
quadruple-precision) then it must emulate them since Intel hardware does
not have them. That would make computations quite slow and you should
consider if you really need that many significant digits.

If GNAT allowed floating-point emulation, you could declare a custom type
of required precision:

   type My_Float is digits 34; -- 112 * lg(2)

Which is advisable anyway in order to make your program portable. Using
built-in types is a bad idea and poor taste in most cases.

Since GNAT does support emulation, at least not with the standard compiler
switches, you need a library for doing this.

-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Fri, 1 May 2015 00:08:36 +0200, Dmitry A. Kazakov wrote:

> Since GNAT does support emulation, at least not with the standard compiler
                           ^^^^ not
> switches, you need a library for doing this.


-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

On 4/30/2015 5:08 PM, Dmitry A. Kazakov wrote:
>
> I don't know what FORTRAN does, but assuming, you are using Intel, the
> longest floating-point available there is not, as you seem assume, of 16
> bytes (quadruple-precision), but so-called extended precision, which is not
> that long (10 bytes, namely).
>

Yes, I do not know either how gfortran does it. I am on intel,
using

>gfortran --version
GNU Fortran (Ubuntu 4.8.2-19ubuntu1) 4.8.2

It might be linked to libquadmath? I see article on this

http://glennklockwood.blogspot.com/2014/02/linux-perf-libquadmath-and-gfortrans.html

"Executive Summary: libquadmath was introduced in GFortran 4.6
which fundamentally changed what the -fdefault-real-8 switch does"

thanks for your input below, will look at it:

> AFAIK, GNAT on x86 maps:
>
> Float to single precision float (4 bytes)
> Long_Float to double precision float (8 bytes)
> Long_Long_Float to extended precision float (10 bytes, effectively)
>
> You can try this program to see how many binary mantissa bits are available
> for each type:
> -----------------------------
> with Text_IO; use Text_IO;
> procedure Test is
> begin
>     Put_Line ("Float mantissa:" & Integer'Image (Float'Machine_Mantissa));
>     Put_Line ("Long Float mantissa:" & Integer'Image
> (Long_Float'Machine_Mantissa));
>     Put_Line ("Long Long Float mantissa:" & Integer'Image
> (Long_Long_Float'Machine_Mantissa));
> end Test;
> --------------------------
> On an Intel machine it will print:
>
> Float mantissa: 24
> Long_Float mantissa: 53
> Long_Long_Float mantissa: 64
>
> If FORTRAN supports longer floating point numbers (e.g.
> quadruple-precision) then it must emulate them since Intel hardware does
> not have them. That would make computations quite slow and you should
> consider if you really need that many significant digits.
>
> If GNAT allowed floating-point emulation, you could declare a custom type
> of required precision:
>
>     type My_Float is digits 34; -- 112 * lg(2)
>
> Which is advisable anyway in order to make your program portable. Using
> built-in types is a bad idea and poor taste in most cases.
>
> Since GNAT does support emulation, at least not with the standard compiler
> switches, you need a library for doing this.
>

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

> On 4/30/2015 5:08 PM, Dmitry A. Kazakov wrote:

> Since GNAT does support emulation, at least not with the standard compiler
> switches, you need a library for doing this.

fyi; from

https://gcc.gnu.org/wiki/GFortran/News

it says:

"GCC now ships with the LGPL-licensed libquadmath library,
which provides for targets with a __float128 datatype
quad-precision mathematical functions. __float128 is
available for targets on 32-bit x86, x86-64 and
Itanium architectures. The libquadmath library is
automatically build on such targets when building the Fortran compiler."

So, I guess the summary of all this is: gfortran supports
this library now directly, but to use it in Ada, more work would
be needed by the user to interface to it (may be gnat needs
also to also add support to allow user code to access this library
as well).

--Nasser

Re: getting same output as gfortran, long_float

[Jacob Sparre Andersen]

Nasser M. Abbasi wrote:

> So, I guess the summary of all this is: gfortran supports this library
> now directly, but to use it in Ada, more work would be needed by the
> user to interface to it (may be gnat needs also to also add support to
> allow user code to access this library as well).

The Ada library PHCpack also contains a quad-double implementation (both
scalar and complax numbers):

   https://github.com/janverschelde/PHCpack/blob/master/src/Ada/Math_Lib/QD/quad_double_numbers.ads

Greetings,

Jacob
-- 
"When we cite authors we cite their demonstrations, not their names"
                                                           -- Pascal

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

On 4/30/2015 5:08 PM, Dmitry A. Kazakov wrote:

> If GNAT allowed floating-point emulation, you could declare a custom type
> of required precision:
>
>     type My_Float is digits 34; -- 112 * lg(2)
>
> Which is advisable anyway in order to make your program portable. Using
> built-in types is a bad idea and poor taste in most cases.
>
> Since GNAT does support emulation, at least not with the standard compiler
> switches, you need a library for doing this.
>

I found that gfortran can do 128 bit floating point without
the use of the compiler switch -fdefault-real-8. Which will
map to similar thing as the above Ada construct:

-------------------
PROGRAM foo
IMPLICIT NONE
REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
PRINT *, x
END PROGRAM
------------------

gfortran -Wall  foo2.f90
./a.out
    1.82898948331047112025871115292829927E-0003

So, I wonder what it would take to integrate this
library to gnat to allow one to do

     type My_Float is digits 34;

and just have it work as is. I assume this requires
support and some changes at the compiler level
(system.max_digits has to change for example)? Now it
is set at 18 and the above gives compiler error.

--Nasser

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Thu, 30 Apr 2015 20:12:47 -0500, Nasser M. Abbasi wrote:

> I found that gfortran can do 128 bit floating point without
> the use of the compiler switch -fdefault-real-8. Which will
> map to similar thing as the above Ada construct:
> 
> -------------------
> PROGRAM foo
> IMPLICIT NONE
> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
> PRINT *, x
> END PROGRAM
> ------------------

I didn't use FORTRAN for decades, but in good old FORTRAN-IV you declare
specific lengths using T*n, e.g.

INTEGER*2
REAL*4
REAL*8

so, I would suggest, maybe naively,

REAL*16

-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

On 5/1/2015 1:52 AM, Dmitry A. Kazakov wrote:
> On Thu, 30 Apr 2015 20:12:47 -0500, Nasser M. Abbasi wrote:
>
>> I found that gfortran can do 128 bit floating point without
>> the use of the compiler switch -fdefault-real-8. Which will
>> map to similar thing as the above Ada construct:
>>
>> -------------------
>> PROGRAM foo
>> IMPLICIT NONE
>> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
>> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
>> PRINT *, x
>> END PROGRAM
>> ------------------
>
> I didn't use FORTRAN for decades, but in good old FORTRAN-IV you declare
> specific lengths using T*n, e.g.
>
> INTEGER*2
> REAL*4
> REAL*8
>
> so, I would suggest, maybe naively,
>
> REAL*16
>

I am no Fortran expert by any means, but yes, REAL*16
will also work. However, modern Fortran suggests to use KIND.

Actually, it is much more complicated. The way I wrote it
above is not the optimal way. One is supposed to call
SELECTED_REAL_KIND(n,e) requesting n significant digits
and e number of digits in the exponent (e is optional).

i.e. one is supposed to write

SELECTED_REAL_KIND(8,3)  specifies real type of
(+-) 0.xxxxxxxx * 10 ^(+-)xxx

If the compiler does not support this, then the compile
will fail. This is in a way similar to Ada's

    type my_type is digits n;

And it is supposed to be portable way of doing things, vs.
using Real*16.  Here is an example:

---------------------
PROGRAM foo
IMPLICIT NONE
INTEGER, PARAMETER :: quad = SELECTED_REAL_KIND(30,10)
REAL(KIND = quad) :: x
x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
PRINT *, x
END PROGRAM
------------------------
>gfortran -Wall  foo2.f90
>./a.out
    1.82898948331047112025871115292829927E-0003

If one would to ask for more digits than the compiler can support,
one gets an -1 from SELECTED_REAL_KIND() call, and this generates
a compile error in the next statement, since KIND=-1 is compile error.

There is much more to this, and it can get really complicated
to understand it all. The bottom line is that using KIND is the
recommended way to define real variables in modern Fortran to obtain
same precision  on different systems.

Here is one link on the subject

http://fortranwiki.org/fortran/show/Real+precision

--Nasser

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Fri, 01 May 2015 02:32:27 -0500, Nasser M. Abbasi wrote:

> Actually, it is much more complicated. The way I wrote it
> above is not the optimal way. One is supposed to call
> SELECTED_REAL_KIND(n,e) requesting n significant digits
> and e number of digits in the exponent (e is optional).

Huh, they finally learned something after half of a century! (:-))

> i.e. one is supposed to write
> 
> SELECTED_REAL_KIND(8,3)  specifies real type of
> (+-) 0.xxxxxxxx * 10 ^(+-)xxx
> 
> If the compiler does not support this, then the compile
> will fail. This is in a way similar to Ada's
> 
>     type my_type is digits n;

Yes, though Ada also mandates that for the precision specified, the
implementation must guarantee certain accuracy of operations. Maybe in the
following 50 years FORTRAN will get that idea as well.

> And it is supposed to be portable way of doing things, vs.
> using Real*16.

Actually REAL*16 is exactly portable. The problem addressed is not
portability, it is rather design being driven by requirements from the
problem space rather than from implementation specific aspects.

-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[robin.vowels]

On Friday, May 1, 2015 at 5:45:23 PM UTC+10, Dmitry A. Kazakov wrote:
> On Fri, 01 May 2015 02:32:27 -0500, Nasser M. Abbasi wrote:
> 
> > Actually, it is much more complicated. The way I wrote it
> > above is not the optimal way. One is supposed to call
> > SELECTED_REAL_KIND(n,e) requesting n significant digits
> > and e number of digits in the exponent (e is optional).
> 
> Huh, they finally learned something after half of a century! (:-))
> 
> > i.e. one is supposed to write
> > 
> > SELECTED_REAL_KIND(8,3)  specifies real type of
> > (+-) 0.xxxxxxxx * 10 ^(+-)xxx
> > 
> > If the compiler does not support this, then the compile
> > will fail. This is in a way similar to Ada's
> > 
> >     type my_type is digits n;
> 
> Yes, though Ada also mandates that for the precision specified, the
> implementation must guarantee certain accuracy of operations. Maybe in the
> following 50 years FORTRAN will get that idea as well.
> 
> > And it is supposed to be portable way of doing things, vs.
> > using Real*16.
> 
> Actually REAL*16 is exactly portable.

No it's not.
Some compilers treat that is an error.
That form is non-standard.

SELECTED_REAL_KIND or a similar modern form is how
precision may be specified in a portable manner.

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Sun, 3 May 2015 17:15:15 -0700 (PDT), robin.vowels@gmail.com wrote:

> On Friday, May 1, 2015 at 5:45:23 PM UTC+10, Dmitry A. Kazakov wrote:

>> Actually REAL*16 is exactly portable.
> 
> No it's not.
> Some compilers treat that is an error.

Portability applies to compilable programs only. You could not compile

   type X is mod 2**64;

either. That does not make it non-portable, IMO.

Portable means: The program exposes same [logically relevant] behavior on
all platforms [of interest].

A program that does not compile on a platform of interest is just an
illegal program. It is neither portable or non-portable. We don't know yet.

> That form is non-standard.

Maybe. Though "FORTRAN" and "standard" sound silly in the same context.

> SELECTED_REAL_KIND or a similar modern form is how
> precision may be specified in a portable manner.

Maybe, but not for the reason of being standard or compilable. The point
was about specifying precision [and accuracy] in the problem space terms
rather than in machine-specific ones.

-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[robin.vowels]

On Monday, May 4, 2015 at 5:21:28 PM UTC+10, Dmitry A. Kazakov wrote:
> On Sun, 3 May 2015 17:15:15 -0700 (PDT), r.nospam@gmail.com wrote:
> 
> > On Friday, May 1, 2015 at 5:45:23 PM UTC+10, Dmitry A. Kazakov wrote:
> 
> >> Actually REAL*16 is exactly portable.
> > 
> > No it's not.
> > Some compilers treat that is an error.
> 
> Portability applies to compilable programs only. You could not compile
> 
>    type X is mod 2**64;
> 
> either. That does not make it non-portable, IMO.
> 
> Portable means: The program exposes same [logically relevant] behavior on
> all platforms [of interest].

Portable means that a program that adheres to the language standard
is capable of running on any system where there is a compiler for that language.

There might be different limits for such things as maximum size of number on
different machines, but nevertheless the program is capable of running on
any particular system provided that it does not exceed such limits.

> A program that does not compile on a platform of interest is just an
> illegal program. It is neither portable or non-portable. We don't know yet.

But in this case, we know that the program is not portable,
because it violates the standard.

> > That form is non-standard.
> 
> Maybe. Though "FORTRAN" and "standard" sound silly in the same context.
> 
> > SELECTED_REAL_KIND or a similar modern form is how
> > precision may be specified in a portable manner.
> 
> Maybe, but not for the reason of being standard or compilable. The point
> was about specifying precision [and accuracy] in the problem space terms
> rather than in machine-specific ones.

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Mon, 4 May 2015 01:53:08 -0700 (PDT), robin.vowels@gmail.com wrote:

> On Monday, May 4, 2015 at 5:21:28 PM UTC+10, Dmitry A. Kazakov wrote:
>> On Sun, 3 May 2015 17:15:15 -0700 (PDT), r.nospam@gmail.com wrote:
>> 
>>> On Friday, May 1, 2015 at 5:45:23 PM UTC+10, Dmitry A. Kazakov wrote:
>> 
>>>> Actually REAL*16 is exactly portable.
>>> 
>>> No it's not.
>>> Some compilers treat that is an error.
>> 
>> Portability applies to compilable programs only. You could not compile
>> 
>>    type X is mod 2**64;
>> 
>> either. That does not make it non-portable, IMO.
>> 
>> Portable means: The program exposes same [logically relevant] behavior on
>> all platforms [of interest].
> 
> Portable means that a program that adheres to the language standard
> is capable of running on any system where there is a compiler for that language.

Consider a system without text output. According to you Ada's Hello World
program is non-portable because such systems exist.

It is a useless definition. For practically any program it would be easy to
present a platform for which it would not compile or work or would expose
some undesired behavior. Which is why it is *not* all behavior, but only
the relevant (contracted) one and not all platforms, but only the intended
ones. There is no such thing as absolute unconditional portability,

> There might be different limits for such things as maximum size of number on
> different machines, but nevertheless the program is capable of running on
> any particular system provided that it does not exceed such limits.

64 bit is just such a limit.

>> A program that does not compile on a platform of interest is just an
>> illegal program. It is neither portable or non-portable. We don't know yet.
> 
> But in this case, we know that the program is not portable,
> because it violates the standard.

A program that violates the standard can still be portable and conversely.

-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[robin.vowels]

On Monday, May 4, 2015 at 8:18:37 PM UTC+10, Dmitry A. Kazakov wrote:
> On Mon, 4 May 2015 01:53:08 -0700 (PDT), r.nospam@gmail.com wrote:
> 
> > On Monday, May 4, 2015 at 5:21:28 PM UTC+10, Dmitry A. Kazakov wrote:
> >> On Sun, 3 May 2015 17:15:15 -0700 (PDT), r.nospam@gmail.com wrote:
> >> 
> >>> On Friday, May 1, 2015 at 5:45:23 PM UTC+10, Dmitry A. Kazakov wrote:
> >> 
> >>>> Actually REAL*16 is exactly portable.
> >>> 
> >>> No it's not.
> >>> Some compilers treat that is an error.
> >> 
> >> Portability applies to compilable programs only. You could not compile
> >> 
> >>    type X is mod 2**64;
> >> 
> >> either. That does not make it non-portable, IMO.
> >> 
> >> Portable means: The program exposes same [logically relevant] behavior on
> >> all platforms [of interest].
> > 
> > Portable means that a program that adheres to the language standard
> > is capable of running on any system where there is a compiler for that language.
> 
> Consider a system without text output. According to you Ada's Hello World
> program is non-portable because such systems exist.

I said no such thing.  You're talking nonsense.

> It is a useless definition. For practically any program it would be easy to
> present a platform for which it would not compile or work or would expose
> some undesired behavior. Which is why it is *not* all behavior, but only
> the relevant (contracted) one and not all platforms, but only the intended
> ones. There is no such thing as absolute unconditional portability,
> 
> > There might be different limits for such things as maximum size of number on
> > different machines, but nevertheless the program is capable of running on
> > any particular system provided that it does not exceed such limits.
> 
> 64 bit is just such a limit.

So?  Did you read what I wrote?
If a program is capable of running on a number of different systems
it is ipso facto portable.

> >> A program that does not compile on a platform of interest is just an
> >> illegal program. It is neither portable or non-portable. We don't know yet.
> > 
> > But in this case, we know that the program is not portable,
> > because it violates the standard.
> 
> A program that violates the standard can still be portable and conversely.

You're talking rubbish.

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Mon, 4 May 2015 06:45:19 -0700 (PDT), robin.vowels@gmail.com wrote:

> On Monday, May 4, 2015 at 8:18:37 PM UTC+10, Dmitry A. Kazakov wrote:
>> On Mon, 4 May 2015 01:53:08 -0700 (PDT), r.nospam@gmail.com wrote:
>> 
>>> On Monday, May 4, 2015 at 5:21:28 PM UTC+10, Dmitry A. Kazakov wrote:
>>>> On Sun, 3 May 2015 17:15:15 -0700 (PDT), r.nospam@gmail.com wrote:
>>>> 
>>>>> On Friday, May 1, 2015 at 5:45:23 PM UTC+10, Dmitry A. Kazakov wrote:
>>>> 
>>>>>> Actually REAL*16 is exactly portable.
>>>>> 
>>>>> No it's not.
>>>>> Some compilers treat that is an error.
>>>> 
>>>> Portability applies to compilable programs only. You could not compile
>>>> 
>>>>    type X is mod 2**64;
>>>> 
>>>> either. That does not make it non-portable, IMO.
>>>> 
>>>> Portable means: The program exposes same [logically relevant] behavior on
>>>> all platforms [of interest].
>>> 
>>> Portable means that a program that adheres to the language standard
>>> is capable of running on any system where there is a compiler for that language.
>> 
>> Consider a system without text output. According to you Ada's Hello World
>> program is non-portable because such systems exist.
> 
> I said no such thing. 

No which thing? Ada RTL without Ada.Text_IO?

>> It is a useless definition. For practically any program it would be easy to
>> present a platform for which it would not compile or work or would expose
>> some undesired behavior. Which is why it is *not* all behavior, but only
>> the relevant (contracted) one and not all platforms, but only the intended
>> ones. There is no such thing as absolute unconditional portability,
>> 
>>> There might be different limits for such things as maximum size of number on
>>> different machines, but nevertheless the program is capable of running on
>>> any particular system provided that it does not exceed such limits.
>> 
>> 64 bit is just such a limit.
> 
> So?  Did you read what I wrote?
> If a program is capable of running on a number of different systems
> it is ipso facto portable.

You said, I quote, "any system where there is a compiler for that
language." So is it "any" or "a number" of systems?

-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[robin.vowels]

On Tuesday, May 5, 2015 at 12:47:17 AM UTC+10, Dmitry A. Kazakov wrote:
> On Mon, 4 May 2015 06:45:19 -0700 (PDT), r.nospam@gmail.com wrote:
> 
> > On Monday, May 4, 2015 at 8:18:37 PM UTC+10, Dmitry A. Kazakov wrote:
> >> On Mon, 4 May 2015 01:53:08 -0700 (PDT), r.nospam@gmail.com wrote:
> >> 
> >>> On Monday, May 4, 2015 at 5:21:28 PM UTC+10, Dmitry A. Kazakov wrote:
> >>>> On Sun, 3 May 2015 17:15:15 -0700 (PDT), r.nospam@gmail.com wrote:
> >>>> 
> >>>>> On Friday, May 1, 2015 at 5:45:23 PM UTC+10, Dmitry A. Kazakov wrote:
> >>>> 
> >>>>>> Actually REAL*16 is exactly portable.
> >>>>> 
> >>>>> No it's not.
> >>>>> Some compilers treat that is an error.
> >>>> 
> >>>> Portability applies to compilable programs only. You could not compile
> >>>> 
> >>>>    type X is mod 2**64;
> >>>> 
> >>>> either. That does not make it non-portable, IMO.
> >>>> 
> >>>> Portable means: The program exposes same [logically relevant] behavior on
> >>>> all platforms [of interest].
> >>> 
> >>> Portable means that a program that adheres to the language standard
> >>> is capable of running on any system where there is a compiler for that language.
> >> 
> >> Consider a system without text output. According to you Ada's Hello World
> >> program is non-portable because such systems exist.
> > 
> > I said no such thing. 
> 
> No which thing? Ada RTL without Ada.Text_IO?

Read the sentence before "I said no such thing.".

Re: getting same output as gfortran, long_float

[Dennis Lee Bieber]

On Fri, 1 May 2015 08:52:43 +0200, "Dmitry A. Kazakov"
<mailbox@dmitry-kazakov.de> declaimed the following:

>On Thu, 30 Apr 2015 20:12:47 -0500, Nasser M. Abbasi wrote:
>
>> I found that gfortran can do 128 bit floating point without
>> the use of the compiler switch -fdefault-real-8. Which will
>> map to similar thing as the above Ada construct:
>> 
>> -------------------
>> PROGRAM foo
>> IMPLICIT NONE
>> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
>> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
>> PRINT *, x
>> END PROGRAM
>> ------------------
>
>I didn't use FORTRAN for decades, but in good old FORTRAN-IV you declare
>specific lengths using T*n, e.g.
>
>INTEGER*2
>REAL*4
>REAL*8
>
	Even that wasn't FORTRAN-IV as I was taught it... F-IV/F66 had

INTEGER
REAL
DOUBLE PRECISION

where INTEGER and REAL used the same size storage.

	The *x notation showed up in extended F-77 versions, as I recall -- and
VAX systems really confused matters with at least four variants (two
different double precision systems; the "normal" one being a single
precision REAL with an additional 32 bits of mantissa [so one could easily
truncate if passing to a single precision function] {F and D float}, the
other with an extended exponent and reduced precision {G float}, and quad
precision {H float})


	The OP is using Fortran-90/95 notation.
-- 
	Wulfraed                 Dennis Lee Bieber         AF6VN
    wlfraed@ix.netcom.com    HTTP://wlfraed.home.netcom.com/

Re: getting same output as gfortran, long_float

[robin.vowels]

On Friday, May 1, 2015 at 4:52:52 PM UTC+10, Dmitry A. Kazakov wrote:
> On Thu, 30 Apr 2015 20:12:47 -0500, Nasser M. Abbasi wrote:
> 
> > I found that gfortran can do 128 bit floating point without
> > the use of the compiler switch -fdefault-real-8. Which will
> > map to similar thing as the above Ada construct:
> > 
> > -------------------
> > PROGRAM foo
> > IMPLICIT NONE
> > REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
> > x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
> > PRINT *, x
> > END PROGRAM
> > ------------------
> 
> I didn't use FORTRAN for decades, but in good old FORTRAN-IV you declare
> specific lengths using T*n, e.g.
> 
> INTEGER*2
> REAL*4
> REAL*8
> 
> so, I would suggest, maybe naively,
> 
> REAL*16

All those statements are non-standard,
and some current compilers treat them as errors.

Re: getting same output as gfortran, long_float

[Dmitry A. Kazakov]

On Thu, 30 Apr 2015 20:12:47 -0500, Nasser M. Abbasi wrote:

> So, I wonder what it would take to integrate this
> library to gnat to allow one to do
> 
>      type My_Float is digits 34;

That would require compiler maintainers.

> and just have it work as is. I assume this requires
> support and some changes at the compiler level
> (system.max_digits has to change for example)? Now it
> is set at 18 and the above gives compiler error.

AFAIK GCC has integrated support for __float128, so if you could live with 

   type IEEE_Float_128 is private;

you would need no changes in the compiler. Arithmetic operations would be
implemented through calls to <quadmath.h> or whatever GCC uses.

-- 
Regards,
Dmitry A. Kazakov
http://www.dmitry-kazakov.de

Re: getting same output as gfortran, long_float

[robin.vowels]

On Friday, May 1, 2015 at 11:12:51 AM UTC+10, Nasser M. Abbasi wrote:
> On 4/30/2015 5:08 PM, Dmitry A. Kazakov wrote:
> 
> > If GNAT allowed floating-point emulation, you could declare a custom type
> > of required precision:
> >
> >     type My_Float is digits 34; -- 112 * lg(2)
> >
> > Which is advisable anyway in order to make your program portable. Using
> > built-in types is a bad idea and poor taste in most cases.
> >
> > Since GNAT does support emulation, at least not with the standard compiler
> > switches, you need a library for doing this.
> >
> 
> I found that gfortran can do 128 bit floating point without
> the use of the compiler switch -fdefault-real-8. Which will
> map to similar thing as the above Ada construct:
> 
> -------------------
> PROGRAM foo
> IMPLICIT NONE
> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad

No it doesn't.  That form is non-standard, and won't compile with some
compilers.
You need to use SELECTED_REAL_KIND or something equivalent.

> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )

And the above constants are double precision, not quadruple precision,
to the numerical result with have no more precision than double.

> PRINT *, x
> END PROGRAM
> ------------------
> 
> gfortran -Wall  foo2.f90
> ./a.out
>     1.82898948331047112025871115292829927E-0003

Re: getting same output as gfortran, long_float

[Jeffrey R. Carter]

On 04/30/2015 02:17 PM, Nasser M. Abbasi wrote:
> 
> with ada.text_io; use ada.text_io;
> procedure foo_2 is
>   x : Long_Long_Float;
> begin
>   x := 12.0 * 0.0001/(1.0 * (1.0 - 0.1)**4 );
>   put_line(Long_Long_Float'image(x));
> end foo_2;

What do you get with this:

With Ada.Text_IO;
procedure Foo_However_Many is
   X : constant := 12.0 * 0.0001 / (1.0 * (1.0 - 0.1) ** 4);
begin -- Foo_However_Many
   Ada.Text_IO.Put_Line (Item => X'Img);
end Foo_However_Many;

? You probably can't get any more significant digits from GNAT than this gives
you. You could also try a type defined as

type Big is digits System.Max_Digits;

-- 
Jeff Carter
"Many times we're given rhymes that are quite unsingable."
Monty Python and the Holy Grail
57

Re: getting same output as gfortran, long_float

[Qun-Ying]

Currently, on my x86_64 Linux with GCC 4.9.2, the maximum digits I could 
define is 18:
type My_float is digits 18;

Far less then what you would expected.

You would need to use some bindings to MPFR/GMP library for more precisions.

Jeffrey R. Carter wrote:
> On 04/30/2015 02:17 PM, Nasser M. Abbasi wrote:
>>
>> with ada.text_io; use ada.text_io;
>> procedure foo_2 is
>>    x : Long_Long_Float;
>> begin
>>    x := 12.0 * 0.0001/(1.0 * (1.0 - 0.1)**4 );
>>    put_line(Long_Long_Float'image(x));
>> end foo_2;
>
> What do you get with this:
>
> With Ada.Text_IO;
> procedure Foo_However_Many is
>     X : constant := 12.0 * 0.0001 / (1.0 * (1.0 - 0.1) ** 4);
> begin -- Foo_However_Many
>     Ada.Text_IO.Put_Line (Item => X'Img);
> end Foo_However_Many;
>
> ? You probably can't get any more significant digits from GNAT than this gives
> you. You could also try a type defined as
>
> type Big is digits System.Max_Digits;
>

Re: getting same output as gfortran, long_float

[Dennis Lee Bieber]

On Thu, 30 Apr 2015 15:27:43 -0700, Qun-Ying <zhu.qunying@gmail.com>
declaimed the following:

>Currently, on my x86_64 Linux with GCC 4.9.2, the maximum digits I could 
>define is 18:
>type My_float is digits 18;
>
>Far less then what you would expected.
>
	More than I'd have expected...

	Classic floating point is considered 7 significant digits for single
precision (4-byte), and 15 significant digits for double precision
(8-byte).
-- 
	Wulfraed                 Dennis Lee Bieber         AF6VN
    wlfraed@ix.netcom.com    HTTP://wlfraed.home.netcom.com/

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

On 4/30/2015 5:12 PM, Jeffrey R. Carter wrote:
> On 04/30/2015 02:17 PM, Nasser M. Abbasi wrote:
>>
>> with ada.text_io; use ada.text_io;
>> procedure foo_2 is
>>    x : Long_Long_Float;
>> begin
>>    x := 12.0 * 0.0001/(1.0 * (1.0 - 0.1)**4 );
>>    put_line(Long_Long_Float'image(x));
>> end foo_2;
>
> What do you get with this:
>
> With Ada.Text_IO;
> procedure Foo_However_Many is
>     X : constant := 12.0 * 0.0001 / (1.0 * (1.0 - 0.1) ** 4);
> begin -- Foo_However_Many
>     Ada.Text_IO.Put_Line (Item => X'Img);
> end Foo_However_Many;
>
> ? You probably can't get any more significant digits from GNAT than this gives
> you. You could also try a type defined as
>
> type Big is digits System.Max_Digits;
>

----------------------------------
gnatmake foo_however_many.adb
gcc-4.8 -c foo_however_many.adb
gnatbind -x foo_however_many.ali
gnatlink foo_however_many.ali
>./foo_however_many
  1.82898948331047096E-03
--------------------------------

And with type Big is digits System.Max_Digits; I get same

----------------------------------
With Ada.Text_IO;
with system;
procedure Foo_However_Many is
    type Big is digits System.Max_Digits;
    X : Big  := 12.0 * 0.0001 / (1.0 * (1.0 - 0.1) ** 4);
begin -- Foo_However_Many
    Ada.Text_IO.Put_Line (Item => X'Img);
end Foo_However_Many;
---------------------------------

./foo_however_many
  1.82898948331047096E-03

thanks,
--Nasser

Re: getting same output as gfortran, long_float

[Jeffrey R. Carter]

On 04/30/2015 03:32 PM, Nasser M. Abbasi wrote:
> 
> With Ada.Text_IO;
> with system;
> procedure Foo_However_Many is
>    type Big is digits System.Max_Digits;
>    X : Big  := 12.0 * 0.0001 / (1.0 * (1.0 - 0.1) ** 4);
> begin -- Foo_However_Many
>    Ada.Text_IO.Put_Line (Item => X'Img);
> end Foo_However_Many;
> 
> ./foo_however_many
>  1.82898948331047096E-03

You could always do

with Ada.Text_IO;
with PragmARC.Rational_Numbers;
with PragmARC.Unbounded_Integers;

procedure Foo_Rational is
   use PragmARC;

   use type Rational_Numbers.Rational;
   use type Unbounded_Integers.Unbounded_Integer;

   One_I  : constant Unbounded_Integers.Unbounded_Integer := +1;
   Twelve : constant Rational_Numbers.Rational := Unbounded_Integers."+" (12) /
One_I;
   P0001  : constant Rational_Numbers.Rational := One_I / (+10_000); -- 0.0001
   P1     : constant Rational_Numbers.Rational := One_I / (+10); -- 0.1

   One : Rational_Numbers.Rational renames Rational_Numbers.One;

   Result : constant Rational_Numbers.Rational := Twelve * P0001 / (One * (One -
P1) ** 4);
begin -- Foo_Rational
   Ada.Text_IO.Put_Line (Item => Rational_Numbers.Image (Result) );
end Foo_Rational;

$ gnatmake -gnatano -gnatwa -O2 -fstack-check -I../RAC foo_rational.adb
gcc-4.6 -c -gnatano -gnatwa -O2 -fstack-check -I../RAC foo_rational.adb
gnatbind -I../RAC -x foo_rational.ali
gnatlink foo_rational.ali -O2 -fstack-check
jrcarter@jcarter-singo-laptop:~/Code$ ./foo_rational
0.0018289894833104709647919524462734339277549154092363968907178783721993598536808413351623228166438042981252857796067672610882487425697302240512117055326931870141746684956561499771376314586191129401005944215820759030635573845450388660265203475080018289894833104709647919524462734339277549154092363968907178783721993598536808413351623228166438042981252857796067672610882487425697302240512117055326931870141746684956561499771376314586191129401005944215820759030635573845450388660265203475080018289894833104709647919524462734339277549154092363968907178783721993598536808413351623228166438042981252857796067672610882487425697302240512117055326931870141746684956561499771376314586191129401005944215820759030635573845450388660265203475080018289894833104709647919524462734339277549154092363968907178783721993598536808413351623228166438042981252857796067672610882487425697302240512117055326931870141746684956561499771376314586191129401005944215820759030635573845450388660265203475080018289894833104709647919

524462

The referenced PragmARC pkgs are available in the beta version of the PragmAda
Reusable Components available from

https://pragmada.x10hosting.com/pragmarc.htm

-- 
Jeff Carter
"Many times we're given rhymes that are quite unsingable."
Monty Python and the Holy Grail
57

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

On 4/30/2015 8:16 PM, Jeffrey R. Carter wrote:

> You could always do
>
> with Ada.Text_IO;
> with PragmARC.Rational_Numbers;
> with PragmARC.Unbounded_Integers;
>
> procedure Foo_Rational is
>     use PragmARC;
>
>     use type Rational_Numbers.Rational;
>     use type Unbounded_Integers.Unbounded_Integer;
>
>     One_I  : constant Unbounded_Integers.Unbounded_Integer := +1;
>     Twelve : constant Rational_Numbers.Rational := Unbounded_Integers."+" (12) /
> One_I;
>     P0001  : constant Rational_Numbers.Rational := One_I / (+10_000); -- 0.0001
>     P1     : constant Rational_Numbers.Rational := One_I / (+10); -- 0.1
>
>     One : Rational_Numbers.Rational renames Rational_Numbers.One;
>
>     Result : constant Rational_Numbers.Rational := Twelve * P0001 / (One * (One -
> P1) ** 4);
> begin -- Foo_Rational
>     Ada.Text_IO.Put_Line (Item => Rational_Numbers.Image (Result) );
> end Foo_Rational;
>
> $ gnatmake -gnatano -gnatwa -O2 -fstack-check -I../RAC foo_rational.adb
> gcc-4.6 -c -gnatano -gnatwa -O2 -fstack-check -I../RAC foo_rational.adb
> gnatbind -I../RAC -x foo_rational.ali
> gnatlink foo_rational.ali -O2 -fstack-check
> jrcarter@jcarter-singo-laptop:~/Code$ ./foo_rational
> 0.0018289894833104709647919524462734339277549154092363968907178783
>72199359853680841335162322816643804298125285779606767261088248742569
>7302240512117055326931870141746684956561499771376314586191129401005
>94421582075903063557384545038866026520347508001828989483310470964791
>95244627343392775491540923639689071787837219935985368084133516232281
>664380429812528577960676726108824874256973022405121170553269318701

> The referenced PragmARC pkgs are available in the beta version of the PragmAda
> Reusable Components available from
>
> https://pragmada.x10hosting.com/pragmarc.htm
>

Thanks, this is very useful. I noticed small difference in output:
side-by-side:

Ada:          0.00182898948331047096479195244627343392775491540..
Mathematica*: 0.00182898948331047096479195244627343392775491540..
gfortran:     0.00182898948331047112025871115292829927

at digit 18, gfortran result is different. But your Ada rational
package gives same result as Mathematica. This tells me your
result is the accurate one !

Mathematica code:
12*1/(10000)/(1*(1 - 1/10)^4) -- do it all in symbolic first
       --->  4/2187

N[%, 100] ; -- ask for 100 digits numerical

humm....I do not know why these are different. I would have expected
gfortran to be accurate for at least 34 digits, not just 18. May be
I am doing something wrong in gfortran. Here is the gfortran again
for reference:

-------------------
PROGRAM foo
IMPLICIT NONE
REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
PRINT *, x
END PROGRAM
------------------

gfortran -Wall  foo2.f90

Re: getting same output as gfortran, long_float

[Jacob Sparre Andersen]

Nasser M. Abbasi wrote:

> Thanks, this is very useful. I noticed small difference in output:
> side-by-side:
>
> Ada:          0.00182898948331047096479195244627343392775491540..
> Mathematica*: 0.00182898948331047096479195244627343392775491540..
> gfortran:     0.00182898948331047112025871115292829927
>
> at digit 18, gfortran result is different. But your Ada rational
> package gives same result as Mathematica. This tells me your
> result is the accurate one !

The problem might be in how quad-double floating point operations work.
(But it could of course also be a plain oldfashioned bug in gfortran.)

Greetings,

Jacob
-- 
»If you're going to have crime,
 it might as well be organized crime.«      -- Lord Vetinari

Re: getting same output as gfortran, long_float

[Waldek Hebisch]

Nasser M. Abbasi <nma@12000.org> wrote:
> 
> Thanks, this is very useful. I noticed small difference in output:
> side-by-side:
> 
> Ada:          0.00182898948331047096479195244627343392775491540..
> Mathematica*: 0.00182898948331047096479195244627343392775491540..
> gfortran:     0.00182898948331047112025871115292829927
> 
> at digit 18, gfortran result is different. But your Ada rational
> package gives same result as Mathematica. This tells me your
> result is the accurate one !

Your previous gfortran result was:

1.82898948331047096479195244627343369E-3

which agrees with others... 

> humm....I do not know why these are different. I would have expected
> gfortran to be accurate for at least 34 digits, not just 18. May be
> I am doing something wrong in gfortran. Here is the gfortran again
> for reference:
> 
> -------------------
> PROGRAM foo
> IMPLICIT NONE
> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
> PRINT *, x
> END PROGRAM

AFAICS you have 34 digit x, but 17 digit constants.  So
it seems that you assign 17 digit number to x, no
wonder there is no gain in accuracy.

-- 
                              Waldek Hebisch
hebisch@antispam.math.uni.wroc.pl 

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

On 5/1/2015 10:39 AM, Waldek Hebisch wrote:
> Nasser M. Abbasi <nma@12000.org> wrote:
>>
>> Thanks, this is very useful. I noticed small difference in output:
>> side-by-side:
>>
>> Ada:          0.00182898948331047096479195244627343392775491540..
>> Mathematica*: 0.00182898948331047096479195244627343392775491540..
>> gfortran:     0.00182898948331047112025871115292829927
>>
>> at digit 18, gfortran result is different. But your Ada rational
>> package gives same result as Mathematica. This tells me your
>> result is the accurate one !
>
> Your previous gfortran result was:
>
> 1.82898948331047096479195244627343369E-3
>
> which agrees with others...
>

Thanks. My mistake. So the command line flag was needed after all!!

>gfortran -Wall -fdefault-real-8  foo2.f90
>./a.out
    1.82898948331047096479195244627343369E-0003

Without the flag, the digit is different

>gfortran -Wall foo2.f90
>./a.out
    1.82898948331047112025871115292829927E-0003

as you say, this is due to D0 in code using double, while
the result "x" is quad.  With the flag in, gfortran will
treate all doubles as "quad", hence the correct result.

I thought with -Wall it will catch all these things, but I think
I need more flags to detect this user error

>> for reference:
>>
>> -------------------
>> PROGRAM foo
>> IMPLICIT NONE
>> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
>> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
>> PRINT *, x
>> END PROGRAM
>
> AFAICS you have 34 digit x, but 17 digit constants.  So
> it seems that you assign 17 digit number to x, no
> wonder there is no gain in accuracy.
>

Yes, needed the command line flag. You are right.

--Nasser

Re: getting same output as gfortran, long_float

[Nasser M. Abbasi]

On 5/1/2015 12:27 PM, Nasser M. Abbasi wrote:
>
> Thanks. My mistake. So the command line flag was needed after all!!
>
>> gfortran -Wall -fdefault-real-8  foo2.f90
>> ./a.out
>      1.82898948331047096479195244627343369E-0003
>
> Without the flag, the digit is different
>
>> gfortran -Wall foo2.f90
>> ./a.out
>      1.82898948331047112025871115292829927E-0003
>
...
> I thought with -Wall it will catch all these things, but I think
> I need more flags to detect this user error
>

Just for completion, there is a flag that checks for this user error:

----------------------------
>gfortran -Wall -Wconversion -Wconversion-extra  foo2.f90
foo2.f90:6.4:

x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
     1
Warning: Conversion from REAL(8) to REAL(16) at (1)
-------------------------------

While this gives clean compile:

gfortran -Wall -Wconversion -Wconversion-extra -fdefault-real-8  foo2.f90

So, Fortran can be made to be "stroger typed", but it I think it takes
more effort than with Ada, which would have detected such
implicit conversions at the language semantics level without
the need of a compiler switches.

--Nasser

Re: getting same output as gfortran, long_float

[robin.vowels]

On Saturday, May 2, 2015 at 3:27:29 AM UTC+10, Nasser M. Abbasi wrote:
> On 5/1/2015 10:39 AM, Waldek Hebisch wrote:
> > Nasser M. Abbasi <nnospam@12000.org> wrote:
> >>
> >> Thanks, this is very useful. I noticed small difference in output:
> >> side-by-side:
> >>
> >> Ada:          0.00182898948331047096479195244627343392775491540..
> >> Mathematica*: 0.00182898948331047096479195244627343392775491540..
> >> gfortran:     0.00182898948331047112025871115292829927
> >>
> >> at digit 18, gfortran result is different. But your Ada rational
> >> package gives same result as Mathematica. This tells me your
> >> result is the accurate one !
> >
> > Your previous gfortran result was:
> >
> > 1.82898948331047096479195244627343369E-3
> >
> > which agrees with others...
> >
> 
> Thanks. My mistake. So the command line flag was needed after all!!
> 
> >gfortran -Wall -fdefault-real-8  foo2.f90
> >./a.out
>     1.82898948331047096479195244627343369E-0003
> 
> Without the flag, the digit is different
> 
> >gfortran -Wall foo2.f90
> >./a.out
>     1.82898948331047112025871115292829927E-0003
> 
> as you say, this is due to D0 in code using double, while
> the result "x" is quad.  With the flag in, gfortran will
> treate all doubles as "quad", hence the correct result.
> 
> I thought with -Wall it will catch all these things, but I think
> I need more flags to detect this user error
> 
> >> for reference:
> >>
> >> -------------------
> >> PROGRAM foo
> >> IMPLICIT NONE
> >> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
> >> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
> >> PRINT *, x
> >> END PROGRAM
> >
> > AFAICS you have 34 digit x, but 17 digit constants.  So
> > it seems that you assign 17 digit number to x, no
> > wonder there is no gain in accuracy.
> >
> 
> Yes, needed the command line flag.

Not really. he correct way is to specify quad precision constants,
e.g.,
0.0001_quad

Re: getting same output as gfortran, long_float

[robin.vowels]

On Friday, May 1, 2015 at 11:40:17 AM UTC+10, Nasser M. Abbasi wrote:
> On 4/30/2015 8:16 PM, Jeffrey R. Carter wrote:
> 
> > You could always do
> >
> > with Ada.Text_IO;
> > with PragmARC.Rational_Numbers;
> > with PragmARC.Unbounded_Integers;
> >
> > procedure Foo_Rational is
> >     use PragmARC;
> >
> >     use type Rational_Numbers.Rational;
> >     use type Unbounded_Integers.Unbounded_Integer;
> >
> >     One_I  : constant Unbounded_Integers.Unbounded_Integer := +1;
> >     Twelve : constant Rational_Numbers.Rational := Unbounded_Integers."+" (12) /
> > One_I;
> >     P0001  : constant Rational_Numbers.Rational := One_I / (+10_000); -- 0.0001
> >     P1     : constant Rational_Numbers.Rational := One_I / (+10); -- 0.1
> >
> >     One : Rational_Numbers.Rational renames Rational_Numbers.One;
> >
> >     Result : constant Rational_Numbers.Rational := Twelve * P0001 / (One * (One -
> > P1) ** 4);
> > begin -- Foo_Rational
> >     Ada.Text_IO.Put_Line (Item => Rational_Numbers.Image (Result) );
> > end Foo_Rational;
> >
> > $ gnatmake -gnatano -gnatwa -O2 -fstack-check -I../RAC foo_rational.adb
> > gcc-4.6 -c -gnatano -gnatwa -O2 -fstack-check -I../RAC foo_rational.adb
> > gnatbind -I../RAC -x foo_rational.ali
> > gnatlink foo_rational.ali -O2 -fstack-check
> > jrcarter@jcarter-singo-laptop:~/Code$ ./foo_rational
> > 0.0018289894833104709647919524462734339277549154092363968907178783
> >72199359853680841335162322816643804298125285779606767261088248742569
> >7302240512117055326931870141746684956561499771376314586191129401005
> >94421582075903063557384545038866026520347508001828989483310470964791
> >95244627343392775491540923639689071787837219935985368084133516232281
> >664380429812528577960676726108824874256973022405121170553269318701
> 
> > The referenced PragmARC pkgs are available in the beta version of the PragmAda
> > Reusable Components available from
> >
> > https://pragmada.x10hosting.com/pragmarc.htm
> >
> 
> Thanks, this is very useful. I noticed small difference in output:
> side-by-side:
> 
> Ada:          0.00182898948331047096479195244627343392775491540..
> Mathematica*: 0.00182898948331047096479195244627343392775491540..
> gfortran:     0.00182898948331047112025871115292829927
> 
> at digit 18, gfortran result is different.

That's because your constants are double precision, not quad precision,
as I said before.

> But your Ada rational
> package gives same result as Mathematica. This tells me your
> result is the accurate one !
> 
> Mathematica code:
> 12*1/(10000)/(1*(1 - 1/10)^4) -- do it all in symbolic first
>        --->  4/2187
> 
> N[%, 100] ; -- ask for 100 digits numerical
> 
> humm....I do not know why these are different. I would have expected
> gfortran to be accurate for at least 34 digits, not just 18.

They cannot be any more accurate than double precision (about 16 digits)
because your constants are double precision.

> May be
> I am doing something wrong in gfortran. Here is the gfortran again
> for reference:
> 
> -------------------
> PROGRAM foo
> IMPLICIT NONE
> REAL(KIND = 16) :: x  !-- kind=16 tells it is double quad
> x = 12.0D0 * 0.0001D0/(1.0D0 * (1.0D0 - 0.1D0)**4 )
> PRINT *, x
> END PROGRAM
> ------------------

Re: getting same output as gfortran, long_float

[Simon Wright]

"Nasser M. Abbasi" <nma@12000.org> writes:

> Now, gfortran has an option called -fdefault-real-8 which I thought
> will map to Ada long_float.

I thought that, with GCC on Intel, Ada Long_Long_Float is the same as
Fortran REAL*10.

But they don't _print_ the same:

A: 1.82898948331047096E-03
F: 1.82898948331047096483E-0003

Re: getting same output as gfortran, long_float

[Georg Bauhaus]

On 30.04.15 23:17, Nasser M. Abbasi wrote:
>
> How to obtain the last result above in Ada?

One option seems to be to use GNAT.SSE as a base, and then overload
operators "*" and the like. This should be good for computations,
not sure about I/O. The package SSE defines a type m128 providing for
SSE operations.

I don't know what the Fortran part of GCC does that the Ada part
couldn't also do, but maybe it's a matter of supply and demand,
supply following demand. As the work seems to have been done
for the Fortran part of GCC, wishful thinking might conclude that
adapting the work to the Ada part's requirements is a reasonably
inexpensive proposition, who knows?

Arbitrary precision might seem another option:
http://web.am.qub.ac.uk/users/j.parker/miscellany/arbitrary/README.arbitrary


The representations chosen by GNAT is seen when passing the switch -gnatRN
to the compiler. For example:

$ gnatmake -gnatwa -gnatR2 foof8.adb
gcc -c -gnatwa -gnatR3 foof8.adb

Representation information for unit Foof8 (body)
------------------------------------------------

for Real'Size use 64;
for Real'Alignment use 8;

for X'Size use 64;
for X'Alignment use 8;
gnatbind -x foof8.ali
gnatlink foof8.ali

Compilation finished at Fri May  1 13:50:13

Where

procedure Foof8
is
    Byte : constant := 8;

    type Real is new Long_Float;
    for Real'Size use 8 * Byte;
    
    --   Dmitry hinted at this:
    --   "If -fdefault-real-8 is given, DOUBLE PRECISION would instead
    --   be promoted to 16 bytes if possible".
    
--   type Double_Precision is digits Natural(112.0 * 0.3)
--   with
--     Size => 16 * Byte;
    
    X : Real;
begin
    X := 0.0;
    if X /= X then
       X := 1.0;
    end if;
end Foof8;