I checked in version 8.6. The code looks to be a bit buggy.
I was able to make your input file work by switching to bilinear elements, manually computing the number of elements/nodes, and fixing an error you had in the input file.
feap **Hill for RVE of alternate layers of 2 Materials **
81 64 2 2 2 4
and
block
cart d d
layer 1 ! layer of material 1 and 2 along direction r (x1 or x)
1 2 1 2 1 2 1 2
quad 4
1 0 0
2 L 0
3 L L
4 0 L
The relevant code is program/blkgen.f and this helps to understand what FEAP is trying to do. Note that nr, ns, nt are then number of nodal divisions that you specified with the CART d d line in your block; so in your case nr = 8 and ns = 8 too.
110 errck = tinput(layer,1,td,5)
if(errck) go to 110
! Layer data
if(pcomp(layer(1),'laye',4)) then
dlayer = nint(td(1))
if(dlayer.eq.1) then
nlay = nr
elseif(dlayer.eq.2) then
nlay = ns
elseif(dlayer.eq.3) then
nlay = nt
endif
errck = palloc(111,'TEMP1',nlay,1)
j = 1
do while(j.le.nlay)
120 errck = tinput(layer,0,td,16)
if(errck) go to 120
do i = j,min(j+15,nlay)
mr(np(111)+i-1) = nint(td(i-j+1))
end do ! i
j = j + 16
end do ! while
go to 110
Note that the direction for LAYEr needs to be {1,2, or 3}. Then the code performs an polling input until you have input data for the number of nodal divisions in the chosen direction, 8 in your case.
The reading happens 16 divisions at a time.
We'll have a look to see what can be done. If you want to accelerate the process, have a look at program/blkgen.f and then program/sblke.f where the material numbers are set using the data that you have entered in the block command. Note that in sblke.f the direction of the layers in in dlayer and the data is in the array ilr. Most likely the error is in the indexing into ilr.