Indeed. Many codes on Netlib share this problem and inspection of the documentation (if any) or the procedure statements to infer the expected subscript range is needed to judge what is (L)
ength supposed to be.
I think the typical example is in the ODE integrators, where you can provide either a scalar absolute tolerance, or a vector. Since these are typically intent(in) you can usually fix this today with an array expression call odeslv(...,[atol],iopt)
. An example is shown here: Has anyone used DVODE and is it good? - #27 by CRquantum
This code was published in 1991, but probably it originated a decade earlier, around 1980 or so. In the accompanying articles, the authors mention running it on an Amdahl 5850. I found a digital copy of an Amdahl FORTRAN 77 Reference Manual, and in Section 3.4 it describes assumed-size semantics. In Section 10.7 it also says,
If an actual argument is an array element name, that element is associated with the initial dummy array element. In this case, the size declared for the dummy array argument must satisfy the formula
size <= (ds+1-as)*(e1)
wheresize
is the storage required by the dummy array argument,ds
is the maximum subscript value that an actual array element can attain,as
is the specified subscript value for the actual array element, ande1
is the storage size occupied by one actual array element.The array in the subroutine must not attempt to access a subscript value that beyond the limit of the array in the invoking program (see chapter 3, “Dummy Array and Actual Array” under paragraph 3.4). […]
According to these rules, the actual argument V(IV)
would amount to a 1-element array. The authors of the code knew this, because in a manual, they say
Note that in the case when there are no ODEs coupled to the PDEs,
V
andVDOT
will be dummy vectors of length one.
but they stop short of saying that the element V(1)
should not be referenced when there are no ODEs (NV = 0
).
Does that mean with Fortran 90, I can correct the code like shown below?
IV = NPDE*NPTS
CALL CSET(NPDE,NPTS,U,WK(I6),WK,WK(I2),WK(I5),NEL,NPTL,WK(I4),
* WK(I12),XBK,IBK,WK(I3),U(IV+1:IV+NV),NV)
This way when NV=0
, the slice U(IV+1:IV+NV)
would be a zero-length array?
EDIT: running a quick check on the program below,
subroutine zerolen(nv,v)
integer nv
double precision v(nv)
print *, nv
v(1) = 3.0
end
C main program
external zerolen
double precision v(6)
C Fortran 90 Solution
nv = 0
call zerolen(nv,v(5+1:5+nv))
end
compiled with -fcheck=bounds
shows the output I wanted to achieve:
At line 5 of file /app/example.f90
Fortran runtime error: Index '1' of dimension 1 of array 'v' above upper bound of 0