Matrix initialization - how to avoid overlapping?

338 views Asked by At

I have a 21x21 matrix which, in majority, holds ones. Certain elements need to be initialized by discrete numbers. The existing solution is a two step initialization: first set all elements to 1 then reassign selected elements to discrete numbers.

This is implemented in a legacy FORTRAN77 code, see below, and it compiles well with gfortran-4.9 (even without warning).

However, I need to compile it on an older server running gfortran-4.4. This fails, together with alternative tools like f2c or the g95 compiler, due to "overlapping initialization".

Now, how can I recode the problem that it avoids this overlap in the special BLOCK DATA program unit? I also think it is not recommended albeit possible with the Intel Fortran Compiler because the order of initialization cannot be guaranteed.

      SUBROUTINE TEST 
      REAL*8 MAT(21,21)
      COMMON /PARAMETERS/ MAT
      END

      BLOCK DATA 
      REAL*8 MAT(21,21)
      COMMON /PARAMETERS/ MAT
      DATA MAT/441*1/
      DATA (MAT(1,J),J=2,19)/
     &     0.971440D0, 0.940444D0, 1, 0.994435D0, 0.708218D0,
     &     0.931484D0, 1.170520D0, 0.990124D0, 1, 1.019530D0,
     &     0.989844D0, 1.002350D0, 0.999248D0, 1.107274D0, 0.880880D0,
     &     0.880973D0, 0.881047D0, 0.881141D0/
      DATA (MAT(2,J),J=3,14)/
     &     1.022740D0, 0.970120D0, 0.945939D0, 0.744954D0, 0.902271D0,
     &     1.084320D0, 1.005710D0, 1.021000D0, 0.944914D0, 0.973384D0,
     &     0.959340D0, 0.945520D0/
      DATA (MAT(3,J),J=4,19)/
     &     0.925053D0, 0.940237D0, 0.849408D0, 0.955052D0, 1.281790D0,
     &     1.5D0, 1, 0.904849D0, 0.897342D0, 0.724255D0,
     &     0.859744D0, 0.855134D0, 0.831229D0, 0.808310D0, 0.784323D0,
     &     0.745171D0/
      DATA (MAT(4,J),J=5,14)/1.022540D0, 0.493148D0, 0.944871D0,
     &     1.144440D0, 3*1, 1.013040D0, 1, 1.00532D0/
      DATA (MAT(5,J),J=8,12)/1.034787D0, 3*1, 1.0049D0/
      DATA (MAT(7,J),J=15,19)/1.008492D0, 1.010124D0, 1.011501D0,
     &     1.012821D0, 1.014089D0/
      DATA (MAT(8,J),J=9,12)/1.1D0, 1, 1.3D0, 1.3D0/
      END
2

There are 2 answers

5
roygvib On BEST ANSWER

If "overlapping initialization" is not allowed, a simple approach may be just to insert n*1 before and after the specific data and fill the remaining part with 1 also, for example,

      BLOCK DATA 
      REAL*8 MAT(21,21)
      COMMON /PARAMETERS/ MAT
      DATA MAT(1,:) / 1,
     &     0.971440D0, 0.940444D0, 1, 0.994435D0, 0.708218D0,
     &     0.931484D0, 1.170520D0, 0.990124D0, 1, 1.019530D0,
     &     0.989844D0, 1.002350D0, 0.999248D0, 1.107274D0, 0.880880D0,
     &     0.880973D0, 0.881047D0, 0.881141D0, 2*1 /
      DATA MAT(2,:) / 2*1,
     &     1.022740D0, 0.970120D0, 0.945939D0, 0.744954D0, 0.902271D0,
     &     1.084320D0, 1.005710D0, 1.021000D0, 0.944914D0, 0.973384D0,
     &     0.959340D0, 0.945520D0, 7*1 /
      DATA MAT(3,:) / 3*1,
     &     0.925053D0, 0.940237D0, 0.849408D0, 0.955052D0, 1.281790D0,
     &     1.5D0, 1, 0.904849D0, 0.897342D0, 0.724255D0,
     &     0.859744D0, 0.855134D0, 0.831229D0, 0.808310D0, 0.784323D0,
     &     0.745171D0, 2*1 /
      DATA MAT(4,:) / 4*1,
     &     1.022540D0, 0.493148D0, 0.944871D0,
     &     1.144440D0, 3*1, 1.013040D0, 1, 1.00532D0, 7*1 /
      DATA MAT(5,:) / 7*1, 1.034787D0, 3*1, 1.0049D0, 9*1 /
      DATA MAT(6,:) / 21*1 /
      DATA MAT(7,:) / 14*1,
     &     1.008492D0, 1.010124D0, 1.011501D0,
     &     1.012821D0, 1.014089D0, 2*1 /
      DATA MAT(8,:) / 8*1, 1.1D0, 1, 1.3D0, 1.3D0, 9*1 /
      DATA MAT(9:21,:) / 273*1 /
      END

As suggested in the comments, another approach is to replace all the DATA statements by array assignments, for example,

      SUBROUTINE init_data
      IMPLICIT NONE
      REAL*8 MAT(21,21), ONE
      COMMON /PARAMETERS/ MAT
      PARAMETER ( ONE = 1.0D0 )
      INTEGER i

      MAT(:,:) = ONE

      MAT(1, 2:19) = [
     &     0.971440D0, 0.940444D0, ONE,        0.994435D0, 0.708218D0,
     &     0.931484D0, 1.170520D0, 0.990124D0, ONE,        1.019530D0,
     &     0.989844D0, 1.002350D0, 0.999248D0, 1.107274D0, 0.880880D0,
     &     0.880973D0, 0.881047D0, 0.881141D0 ]
      MAT(2, 3:14) = [
     &     1.022740D0, 0.970120D0, 0.945939D0, 0.744954D0, 0.902271D0,
     &     1.084320D0, 1.005710D0, 1.021000D0, 0.944914D0, 0.973384D0,
     &     0.959340D0, 0.945520D0 ]
      MAT(3, 4:19) = [
     &     0.925053D0, 0.940237D0, 0.849408D0, 0.955052D0, 1.281790D0,
     &     1.5D0,      ONE,        0.904849D0, 0.897342D0, 0.724255D0,
     &     0.859744D0, 0.855134D0, 0.831229D0, 0.808310D0, 0.784323D0,
     &     0.745171D0 ]
      MAT(4, 5:14) = [ 1.022540D0, 0.493148D0, 0.944871D0,
     &     1.144440D0, (ONE,i=1,3), 1.013040D0, ONE, 1.00532D0 ]
      MAT(5, 8:12) = [ 1.034787D0, (ONE,i=1,3), 1.0049D0 ]
      MAT(7, 15:19) = [ 1.008492D0, 1.010124D0, 1.011501D0,
     &                  1.012821D0, 1.014089D0 ]
      MAT(8, 9:12) = [ 1.1D0, ONE, 1.3D0, 1.3D0 ]
      END

which seems somewhat more neat to me. Though I have confirmed that this code works for gfortran 4.4.7, it may be necessary to use (/.../) instead of [...] for more old compilers, such that (/ 1.0d0, 2.0d0, ... /). Please also note that in the second approach we need to call init_data at some point of the program, in contrast to BLOCK DATA where the data are automatically initialized upon execution of the program.

0
I. Lagaris On

One may construct a routine that will be called just once. The routine will read the data for the matrix from a prepared file. That will certainly solve your problem with any compiler.