Error #6901: The decimal constant was too large when converting to an integer, and overflow occurred

Dear all

A quick question. I am writing a Fortran subroutine which contains some large elements in p array as below,

  subroutine expmvtay1() 
  use iso_fortran_env
  implicit none
  integer(int32), parameter, dimension(3) :: p=[ &
  8320987112741391580056396102959641077457945541076708813599085350531187384917164032 &
, 507580213877224835833540161373088490724281389843871724559898414118829028410677788672 &
, 31469973260387939390320343330721249710233204778005956144519390914718240063804258910208]
  end subroutine expmvtay1

However, it compile with error messages below (Intel Fortran):

 error #6901: The decimal constant was too large when converting to an integer, and overflow occurred.   [8320987112741391580056396102959641077457945]

Change int32 to int64 or int128 does not seem to work.
So how do I make the above Fortran subroutine compile and work correctly?

Many thanks in advance!

PS1.
I can think of just make this array p as real type instead of integer. Then the below code can compile without problem

  subroutine expmvtay1() 
  use iso_fortran_env
  implicit none
  real(real64), parameter, dimension(3) :: p=[ & 
8320987112741391580056396102959641077457945541076708813599085350531187384917164032.0_real64 &
, 507580213877224835833540161373088490724281389843871724559898414118829028410677788672.0_real64 &
, 31469973260387939390320343330721249710233204778005956144519390914718240063804258910208.0_real64]
  end subroutine expmvtay1

PS2.
The problem arise from that I am translating a matrix exponential solver in this paper

written in matlab as below
https://personales.upv.es/joalab/software/expmvtay1.m
into Fortran.

The Matlab code is very simple, less than 60 line, however, in the Matlab code, there is a 63-element integer (or real) array p with some very big elements, like below

p=[1, 1, 2, 6, 24, 120, 720, 5040, 40320, 362880, 3628800, 39916800, 479001600, 6227020800, 87178291200, 1307674368000, 20922789888000, 355687428096000, 6402373705728000, 121645100408832000, 2432902008176640000, 51090942171709440000, 1124000727777607680000, 25852016738884978212864, 620448401733239409999872, 15511210043330986055303168, 403291461126605650322784256, 10888869450418351940239884288, 304888344611713836734530715648, 8841761993739700772720181510144, 265252859812191032188804700045312, 8222838654177922430198509928972288, 263130836933693517766352317727113216, 8683317618811885938715673895318323200, 295232799039604119555149671006000381952, 10333147966386144222209170348167175077888, 371993326789901177492420297158468206329856, 13763753091226343102992036262845720547033088, 523022617466601037913697377988137380787257344, 20397882081197441587828472941238084160318341120, 815915283247897683795548521301193790359984930816, 33452526613163802763987613764361857922667238129664, 1405006117752879788779635797590784832178972610527232, 60415263063373834074440829285578945930237590418489344, 2658271574788448529134213028096241889243150262529425408, 119622220865480188574992723157469373503186265858579103744, 5502622159812088456668950435842974564586819473162983440384, 258623241511168177673491006652997026552325199826237836492800, 12413915592536072528327568319343857274511609591659416151654400, 608281864034267522488601608116731623168777542102418391010639872, 30414093201713375576366966406747986832057064836514787179557289984, 1551118753287382189470754582685817365323346291853046617899802820608, 80658175170943876845634591553351679477960544579306048386139594686464, 4274883284060025484791254765342395718256495012315011061486797910441984, 230843697339241379243718839060267085502544784965628964557765331531071488, 12696403353658276446882823840816011312245221598828319560272916152712167424, 710998587804863481025438135085696633485732409534385895375283304551881375744, 40526919504877220527556156789809444757511993541235911846782577699372834750464, 2350561331282878906297796280456247634956966273955390268712005058924708557225984, 138683118545689864933221185143853352853533809868133429504739525869642019130834944, 8320987112741391580056396102959641077457945541076708813599085350531187384917164032, 507580213877224835833540161373088490724281389843871724559898414118829028410677788672, 31469973260387939390320343330721249710233204778005956144519390914718240063804258910208];

That Fortran subroutine example is a minimal example which try to contain this p array.

Looks like a list of factorials of 0 to 62. If you actually need very high precision you would have to look at some of the arbitrary precision libraries available in Fortran; but I suspect (I did not look at the main code yet) that you can just use the gamma(3f) function as show to generate approximations to the factorial values as float values.

A little sample program using brute-force products and comparing them to the results of the gamma function as well as showing the values as integers up to huge(0,kind=int32|int64) is attached if you want to verify what the series is and look at what accuracy results with commonly available floats. Note real128 is not shown but could be (nvfortran does not support real128; there might be others as well that do not).

      program demo_gamma
      use, intrinsic :: iso_fortran_env, only : wp=>real64, int64
      implicit none
      real :: x, xa(4)
      integer :: i, j

        ! gamma() is related to the factorial function
        do i = 1, 171
           ! check value is not too big for default integer type
           if (factorial(i)  <=  huge(0)) then
              write(*,*) i, nint(factorial(i)), 'integer'
           elseif (factorial(i)  <=  huge(0_int64)) then
              write(*,*) i, nint(factorial(i),kind=int64),'integer(kind=int64)'
           else
              write(*,*) i, factorial(i) , 'user factorial function'
              write(*,*) i, product([(real(j, kind=wp), j=1, i)]), 'product'
              write(*,*) i, gamma(real(i + 1, kind=wp)), 'gamma directly'
           endif
        enddo

      contains
      function factorial(i) result(f)
      !  GAMMA(X) computes Gamma of X. For positive whole number values of N the
      !  Gamma function can be used to calculate factorials, as (N-1)! ==
      !  GAMMA(REAL(N)). That is
      !
      !      n! == gamma(real(n+1))
      !
      integer, intent(in) :: i
      real(kind=wp) :: f
        if (i  <=  0) then
           write(*,'(*(g0))') '<ERROR> gamma(3) function value ', i, ' <= 0'
           stop '<STOP> bad value in gamma function'
        endif
        f = anint(gamma(real(i + 1,kind=wp)))
      end function factorial

      end program demo_gamma
1 Like

The shorter answer is the gamma(3f) function is the Fortran intrinsic that can generate factorial approximations and it is elemental, so to precompute the factorials (probably to improve performance rather than calculating them every time) you can define the vector P with the following

      program demo_gamma
      use, intrinsic :: iso_fortran_env, only : wp=>real64
      implicit none
      integer :: i
      real(kind=wp),parameter :: p(63) = gamma([(real(i + 1,kind=wp),i=0,62)])
      write(*,*)p
      end program demo_gamma
2 Likes

Thank you so much @urbanjost ! So awesome! :100:

I have another small question. In that Matlab code there are things like

s^m

where both s and m are integers. For example s=6, m=63, so s^m=6^63=10^49 or so.
Thing is Matlab have no problem calculating 6^63,

Now, in Fortran, here both s and m are not big intergers, they are at most 60, so I may casually define them as integer(4). However in this example, the result s^m=6^63 cannot be represented by integer(4), integer(8) perhaps even higher should be used. So if I define s and m as integer(4) I got wrong results for s^m.

How would you handle case like that?
I mean, do you manually change

s**m

to

int(s,kind=8)**m

or

real(s,kind=8)**m

if you suspect that s**m may overflow?

Many thanks!

PS.
This is perhaps somehow similar with some other posts before

This is probably a good example of how to use the power of Fortran kinds to specify how to select a kind that needs a specified range of values, but last I knew Matlab converts every numeric value to a doubleprecision multi-dimensional array unless you explicitly use the arbitrary precision routine(s) like VPA() or something like that.

So if 60^60 really is your upper bound that you need to be able to process that is a 107 digit value; which is too big for 32 big integers and floats, and too big for 64 bit integers but will
fit in a 64-bit float. The HUGE() function is an easy way to get the maximum value a type will hold. So changing s^m to

real(s,kind=wp)**m

would work. And make sure the resulting value is returned into a type of sufficient size to hold the value without overflow as well.

program demo_overflow
!60^60=> 48873677980689257489322752273774603865660850176000000000000000000000000000000000000000000000000000000000000
! 107 digits
use, intrinsic :: iso_fortran_env, only : wp=>real64
!implicit real(kind=wp)(a-z)
implicit none
character(len=:),allocatable :: value
real :: float
real(kind=wp) :: bigfloat

   write(*,*)'biggest float=',huge(float)
   write(*,*)'biggest bigfloat=',huge(bigfloat)
   value='48873677980689257489322752273774603865660850176000000000000000000000000000000000000000000000000000000000000'
   read(value,*)float      ! overflow
   write(*,*)float         ! error

   read(value,*)bigfloat
   write(*,*)bigfloat

   ! assuming worse case really is 60**60
   write(*,*)60**60        ! error; result is much bigger than largest default integer
   !write(*,*)60_wp**60    ! error (used real kind on an integer, no decimal)
   write(*,*)60.0_wp**60
end program demo_overflow
 biggest float=   3.40282347E+38
 biggest bigfloat=   1.7976931348623157E+308
         Infinity
   4.8873677980689261E+106
 -2147483648
   4.8873677980689261E+106
1 Like

The LOG_GAMMA intrinsic is useful when you need to calculate the ratio of gamma functions or factorials.

You can calculate GAMMA(1000.0)/GAMMA(999.0) as EXP(LOG_GAMMA(1000.0)-LOG_GAMMA(999.0)). There may be some loss of precision, but it is much less likely to overflow.

$ cat log_gamma_01.f90
  write(*,*) 'Single precision: ',EXP(LOG_GAMMA(1000.0)-LOG_GAMMA(999.0))
  write(*,*) 'Double precision: ',EXP(LOG_GAMMA(1000.0D0)-LOG_GAMMA(999.0D0))
  end

$ gfortran -o log_gamma_01.exe log_gamma_01.f90

$ ./log_gamma_01.exe
Single precision:    998.983521
Double precision:    998.99999999951149






2 Likes

In MATLAB every numeric value is double by default, so keep that in mind when porting your code to Fortran.

1 Like