A weird problem when using vdsqr

[junuylia]

Hi,

I've got a very weird problem. I used the vdsqr function in a loop with k as index, the program was like this:

do k = 1, kk


call vdSQR(32, XI(K,:), RSQ(K, :))

enddo

kk is a very large number, xi and rsq are two matrix with kk rows and 32 columns.

The problem comes when k reaches a certain number, say 10189. At this iteration, after the call, the xi will change its value, and k will become 0. The other iterations when k is less than this number seem to be right.

Also, if I use l instead of k, say,

do l = 1, kk

(my program)
enddo

the problem doesn't show up.

I'm just curious is there any reported problem like this, or is it just my programing mistakes? And any hint on what happend? Thank you very much!

June

[mecej4]

You are using a noncontiguous array section where an array argument is expected.

When you do so, the actual argument is regarded as an expression and the corresponding dummy argument in the called subroutine should not be defined or redefined, and must not have intent out or inout (see Sec. 6.13, Metcalf, Reid and Cohen, Fortran 95/2003 Explained, OUP, 2004).

However, for the routine that you called, the third argument is intent out. Thus, there is a problem, and it can be fixed by using an intermediate 1-D array to hold the third argument. The compiler may, depending on the options in effect, create an array temporary, use the temporary as the argument, and copy the returned values back to the array section in proper sequence.

In other words,

[fortran]real :: sqxi(32)
...
do k = 1, kk

   call vsSQR(32, XI(K,:), sqxi)
   RSQ(K, :) = sqxi

enddo[/fortran]

To see what is happening with your program, I think we need to see the actual code.

[junuylia]

Thanks for the reply, but I don't think this is the issue, or the previous iterations would not go. I did some testing and believe that it has something to do with the memory. Since the array is too large, say 17073*32 = 546336 numbers with type real(8), the system might have trouble allocating such a large array. So when it reaches the limit, the calculation will be crazy.

I think the solution would be that I have to try to avoid allocating this array. If it has to be used, I would try to divide it in parts.


Also, since this code is included in a subroutine, when I defined the array, I used

real(8), dimension(sum(ll)) :: rsq

where ll is a parameter. But this seems to be not appropriate. When I replace sum(ll) with a input argument, the programs runs fine.

[TimP]

As mecej4 said, we would need actual code and specifics of your system to do more than guess. The default Intel OpenMP thread stack size is 2MB on 32-bit systems, 4MB on 64-bit, so you may want to check that if you allocate thread private data.