 # Forall loop confusion

The following function describes a map between two vectors of dimension n to a vector of dimension n^2. Here, vector is a type containing a rank-1 pointer, “values”, and a logical scalar, “tmp”.

``````module c_algebra

type vector
private
complex(dp), dimension(:), pointer :: values => null()
logical                            :: tmp    =  .false.
end type vector

function vec_dir_vec(v1,v2 ) result(v3)
implicit none
type(vector), intent(in)  :: v1,v2
type(vector), intent(out) :: v3
integer                   :: n1, n2, i, j

n1=size(v1); n2=size(v2)
allocate(v3%values(n1 * n2))

!    forall ( i = 1:n1, j = 1:n2 )
do i=1, n1
do j=1, n2
v3%values(2*(i-1)+j) = v1%values(i) * v2%values(j)
end do
end do
!    end forall

v3%tmp=.true.
call cleanup( v1 , .true. )
call cleanup( v2 , .true. )
end function vec_dir_vec
end module c_algebra
``````

Compiling and running with portland 6.1-5 on a x86-64 machine using RHEL produces the following results for v1%values=(/cmplx(1,0),cmplx(2,0/) and v2%values=(/(3,0),(4,0)/). With a forall loop, v3%values=((2,0),(4,0),(2,0),(4,0)), whereas with the nested do loops, the correct result of ((1,0),(2,0),(2,0),(4,0)) is produced.

Now I’m pretty sure the two should be identical. Can anybody help me here?

Further, should the following work, where a1 and a2 are integer arrays, larray is a logical array (all size n), and flag is a logical scalar?

``````larray=a1.eq.a2
flag=.true.
forall ( i=1:n )
if (.not.larray(i)) then
flag=.false.
return
end if
end forall
``````

Or is there perhaps an inbuilt function to discover whether arrays are identical?

``````forall( i = 1:n1, j = 1:n2)
v3(i, j) = expression
end forall
``````

mapping to:

``````do i = 1:n1
do j = 1:n2
v3(i,j) = expression
end do
end do
``````

All right hand side expressions get evaluated before the assignments, and may be calculated out of order or in parallel. The left hand side array indices are also calculated in any order, and must not overlap. Many-one assignments (as your code has) are not allowed.

The bottom line is that you need to use the do loops here.

Hope this helps,
Mat

Just to clear up any confusion. There was a typo here. The forall loop was supposed to be

``````forall ( i = 1:n1; j = 1:n2 )
v3%values(n2*(i-1)+j) = v1%values(i) * v2%values(j)
end forall
``````

… which produces the same spurious answer.

As you can see, this is actually one-to-one, in the sense that any pair of indices (i,j) in (1:n1, 1:n2) will produce a unique index on the LHS. Whilst I can understand that it may be difficult for the compiler to distinguish this from a many-one assignment, it is not such, so I think the standard rather expects it to work. Certainly had I not tested this carefully, I’d not have spotted it. At the very least, should not the compiler flag such things? I’m pretty sure the code is legal, and it’s a problem if it simply doesn’t work.

Hi Che,

Thanks for the clarification. I’m able to reproduce the error now and agree that it’s compiler problem. Note that it seems to only occur when a derived type’s allocatable array is used on the left-hand side of the assignment. I have filed a technical problem report (TPR#4045) and will have our compiler team investigate.

• Mat

Many thanks indeed.