Re: Possible explaination why large_real_kind_form_io_2.f90 is failing on PPC/OSX

[Tim Prince <tprince at myrealbox dot com>]

Dominique Dhumieres wrote:
> > What is the output of -fdump-tree-original?
> >     
>
> MAIN__ ()
> {
>   _gfortran_set_std (70, 127, 0);
>   {
>     struct __st_parameter_dt dt_parm.0;
>
>     dt_parm.0.common.filename = "print_nearest_16.f90";
>     dt_parm.0.common.line = 1;
>     dt_parm.0.common.unit = 6;
>     dt_parm.0.format = "(1pg60.40)";
>     dt_parm.0.format_len = 10;
>     dt_parm.0.common.flags = 4096;
>     _gfortran_st_write (&dt_parm.0);
>     {
>       static real16 C.948 = 1.0e+0;
>
>       _gfortran_transfer_real (&dt_parm.0, &C.948, 16);
>     }
>     _gfortran_st_write_done (&dt_parm.0);
>   }
>
> This seems ver suspicious. I have changed the code to:
>
> print *, precision(1.0_16)
> print '(1pg60.40)', nearest(1.0_16,1.0_16)
> print *, nearest(1.0_16,1.0_16)-1.0_16
> end
>
> The output is:
>
>           31
>               1.0000000000000000000000000000000000000000    
>   2.4651903288156618919116517665087070E-0032
>
> and the dump:
>
> MAIN__ ()
> {
>   _gfortran_set_std (70, 127, 0);
>   {
>     struct __st_parameter_dt dt_parm.0;
>
>     dt_parm.0.common.filename = "print_nearest_16.f90";
>     dt_parm.0.common.line = 1;
>     dt_parm.0.common.unit = 6;
>     dt_parm.0.common.flags = 128;
>     _gfortran_st_write (&dt_parm.0);
>     {
>       static int4 C.948 = 31;
>
>       _gfortran_transfer_integer (&dt_parm.0, &C.948, 4);
>     }
>     _gfortran_st_write_done (&dt_parm.0);
>   }
>   {
>     struct __st_parameter_dt dt_parm.1;
>
>     dt_parm.1.common.filename = "print_nearest_16.f90";
>     dt_parm.1.common.line = 2;
>     dt_parm.1.common.unit = 6;
>     dt_parm.1.format = "(1pg60.40)";
>     dt_parm.1.format_len = 10;
>     dt_parm.1.common.flags = 4096;
>     _gfortran_st_write (&dt_parm.1);
>     {
>       static real16 C.950 = 1.0e+0;
>
>       _gfortran_transfer_real (&dt_parm.1, &C.950, 16);
>     }
>     _gfortran_st_write_done (&dt_parm.1);
>   }
>   {
>     struct __st_parameter_dt dt_parm.2;
>
>     dt_parm.2.common.filename = "print_nearest_16.f90";
>     dt_parm.2.common.line = 3;
>     dt_parm.2.common.unit = 6;
>     dt_parm.2.common.flags = 128;
>     _gfortran_st_write (&dt_parm.2);
>     {
>       static real16 C.952 = 2.46519032881566189191165176650870696772877010972e-32;
>
>       _gfortran_transfer_real (&dt_parm.2, &C.952, 16);
>     }
>     _gfortran_st_write_done (&dt_parm.2);
>   }
>
> If I understand correctly:
>
>   
> > Here, C.991 is a value from MPFR.  D.992 is going to come from GCC.
> >     
>
> It seems the problem is coming from the gfortran/mpfr side.
>
> Dominique
>
>   
Who is to say what the correct relationship is between 
non-IEEE-compliant IBM style real(16) (107 bits maximum precision) and 
IEEE-854 style 107 bits precision, if that is what mpfr uses?
In case it is of interest, here is the ifort 9.1 emulation of IEEE real(16):
         33
              1.000000000000000000000000000000000000000
 1.925929944387235853055977942584927E-0034

ifort does actually call its library function for nearest(); evidently, 
even the effort to emulate IEEE arithmetic is not entirely successful.