Extended precision

Sorry for this question but unfortunately I cannot find any 'clear' information... Is there any ability to work with extended precision (real*10)?

I wrote the following: real*10 A

then I try to use sqrt for A.

But unfortunately it's not possible.

Could anybody give me an advice?:)

Thanks!

The following simple code works ok for me real10 a,b ! a = 8 b = sqrt (a) write (,*) a,b, b/2 end

Quoted from JohnCampbell The following simple code works ok for me real10 a,b ! a = 8 b = sqrt (a) write (,*) a,b, b/2 end

Hm, thanks) will try...

I used 'Dsqrt'....:)

The following code might better demonstrate the precision of the alternative real kinds available in ftn95 real10 a,b real8 c,d real4 e,f ! write (,*) '! sqrt_two = 1.4142135623730950488016887242096980_16 = SQRT (2)'

a = 8 
b = sqrt (a) 
write (*,2010) '*10',a,b, b/2, (b/2)**2
! 
c = 8 
d = sqrt (c) 
write (*,2008) '*8 ',c,d, d/2, (d/2)**2 
! 
e = 8 
f = sqrt (e) 
write (*,2004) '*4 ',e,f, f/2, (f/2)**2 
2010 format (a4, f27.19,3f27.19)
2008 format (a4, f24.16,3f27.16)
2004 format (a4, f16.8, 3f27.8)
end

I've pushed the output precision until sqrt(2) is wrong in the last digit.

Note that b/2 retains the real*10 precision, as integer 2 is an exact binary representation of 2.0_4, while 0.1_2 is different from 0.1_3 or 0.1_4.

real10 does offer better precision than real8, but I have found it difficult to maintain the precision, especially with mixed kind calculations. If your calculation is not stable with real8, then going to real10 might not be the best solution. Changing the calculation formulation, or even the order of calculation might improve round-off growth. For infinite series with alternating sign, it can be useful to do the calculation in pairs to overcome subtracting similar values, although the round-off improvement is often marginal.

John

I had the code waiting for extended precision for ages and i gave it a try today.

Just for info, the machine ESP for real10 is in the region of 510^-20 which can be taken from this simple code

   real*10 Eps, tol 
   EPS=1.D0
10 Eps=Eps/2.D0
   Tol=1.d0+Eps
   if (Tol.GT.1.d0) goto 10
   Write (*,*)' Eps REAL*10=', Eps
   end

but what interesting for me is that the exponent can go beyond the crazy values like 10^4000:

   a=1
   do i=1,20000
   a=a*2
   write (*,*)'Max=', a
   enddo

Question is though: how to use this extended range of exponents in the codes because compiler swears with the message '*** Floating point number out of range *** Compilation abandoned you moron' for that kind of lines

if(SumP.GT.1.d3000) goto  2

Use those annoying REAL_KINDs 3? integer, parameter :: rkind=selected_real_kind(18,4931) How correctly set it up?

If you declare REAL(3)::x then the compiler will provide a useful comment to tell you what to use in place of 3.

Thanks. I also found something relevant in the Help, but will this fix the error i've got?

I don't think so. You will need to write code that avoids overflow. You can probably trap the exception if it occurs at run time but I don't see how that will help.

But that is not an error and not an under/overflow! This is completely legitimate code which got these kind of too large or too small numbers due to tricky method of solution! And these sort of crazy numbers were in the FP coprocessors since IBM PC XT, we just never used them because compilers did not allow that

Maybe you should use constants with the corresponding KIND.

e.g. x = 1.0_3

Alas, compiler still gives an error that floating point number out of rangeBad value for exponent. I suspect the limit in compiler is done for real*8 numbers and is explicitly wired somewhere in the source

Dan,

I'd like to see the example which fails with bad exponent value, including declarations. I think you have to be careful with constant definition, to ensure the implied precision. The following example code fails for the last two lines, but the previous 3 do work.

module kind_precision
   integer, parameter :: int_1byte   = selected_int_kind (2)         ! 1
   integer, parameter :: int_2byte   = selected_int_kind (4)         ! 2
   integer, parameter :: int_4byte   = selected_int_kind (9)         ! 3
   integer, parameter :: int_8byte   = selected_int_kind (18)        ! 4
   integer, parameter :: real_float  = selected_real_kind (6,37)     ! 1
   integer, parameter :: real_double = selected_real_kind (15,307)   ! 2
   integer, parameter :: real_long   = selected_real_kind (18,4931)  ! 3
end module kind_precision

use kind_precision
!
   integer ( kind=int_4byte ) :: i,k, last_k
   character kind_name*10
   real ( kind=real_double ) :: x
   real ( kind=real_long )   :: y
!
!  confirm integer precision
      write (*,*) ' '
      write (*,*) 'Test of integer KIND'
      last_k = -99
      do i = 50,-50,-1
         k = selected_int_kind (i)
         if (k == last_k) cycle
         kind_name = ' undefined'
         if ( k == int_1byte ) kind_name = '1 byte'
         if ( k == int_2byte ) kind_name = '2 byte'
         if ( k == int_4byte ) kind_name = '4 byte'
         if ( k == int_8byte ) kind_name = '8 byte'
         write (*,*) kind_name, ' precision =', i, ' kind =',k
         last_k = k
      end do
!
!  confirm real precision
      write (*,*) ' '
      write (*,*) 'Test of Real precision KIND'
      last_k = -99
      do i = 50,-50,-1
         k = selected_real_kind (i,1)
         if (k == last_k) cycle
         kind_name = ' undefined'
         if ( k == real_float  ) kind_name = '4 byte'
         if ( k == real_double ) kind_name = '8 byte'
         if ( k == real_long   ) kind_name = '10 byte'
         write (*,*) kind_name, ' precision =', i, ' kind =',k
         last_k = k
      end do
!
!  confirm real exponent
      write (*,*) ' '
      write (*,*) 'Test of Real exponent KIND'
      last_k = -99
      do i = 5000,-5000,-1
         k = selected_real_kind (1,i)
         if (k == last_k) cycle
         kind_name = ' undefined'
         if ( k == real_float  ) kind_name = '4 byte'
         if ( k == real_double ) kind_name = '8 byte'
         if ( k == real_long   ) kind_name = '10 byte'
         write (*,*) kind_name, ' exponent =', i, ' kind =',k
         last_k = k
      end do
!
!   test real constants
!
      write (*,*) ' '
      x = 1.0e300_real_double   ; write (*,*) 'x=',x
      y = 1.0e3000_real_long    ; write (*,*) 'y=',y
      x = 1.0d300               ; write (*,*) 'x=',x
!      x = 1.0e300               ; write (*,*) 'x=',x
!      y = 1.0e3000              ; write (*,*) 'y=',y
!
      end

You might like to consider defining your own Fortran 2008 ISO_C_BINDING module, as this provides some standardisation to naming kinds. We have to suffer because Fortran 90 did not recommend that KIND = number_of_bytes. You can try to be too smart !! I'd like to know the Fortran 90+ implementations where this could not have been done. Were there any 60-bit Cyber computers still around in 1990 ?

John

x = 1d400 fails but x = 1e400_3 is OK.

Thanks John and Paul. So

if(aa.gt.2e-4000_3) write(,) 'something'

works but

if(aa.gt.2d-4000_3) write(,) 'something'

DOES NOT work!? THAT....that...that's TOTALLY counterintuitive!

And inconvenient. When you program you adjust the value of constant and use E till 10^{37, then you use D till 10}309 and then again E??? Inconvenient, confusing, counterintuitive. That will end up with the loss of the whole day for any programmer in future. Or more, because i have got this forum and the pros here, others may lose much more than that. Is this how Standard precludes to do that? If these letter serve as a switch for different implementations of the intrinsic functions when you change precision then for extended precision Standard should allow all three different letters E (for 'exponent'), D(from 'Double precision') and, say, X for 'extended precision'! Or at least compiler should give the diagnostics and recommendation (giving recommendation would be next good thing for future Fortran compilers)

I guess that the D is for double precision so in FTN95 D is equivalent to E with _2 as the KIND.