Fortran frozen after iterations in a loop

[pedrohcgs]

Hi,

The code I have attached is freezen at the 204th loop. 

I have no idea the reason for that. 

When I set the parameter nr=100, everything works fine, and when I set nr=500, the code freezes at the 1275th loop.

Since I am new at this, could anyone provide me any help?

I am using Intel(R) Visual Fortran Intel(R) 64 Compiler XE on Intel(R) 64, version 12.1.2 , with IMSL 6.0

Thank you very much!

Pedro

[bmchenry]

since you have a lot going on in the code, it is asking a wee bit too much to have the group walk through it to help you.

what is happening at the 1275th loop?

i'd suggest you put a breakpoint in to stop at the 1274th loop and then walk through the loop.

I would bet it will become evident where the freeze is occuring and then you can correct the issue.

[Steven_L_Intel1]

I can see that RNPOI is getting into a tight loop at that call. It will take me some time to dig into this further, but I will do so.

One thing I see is that you call the D_xxxxx names, such as D_RNNOA. While this is not wrong, it is not the way you should be using these routines.  Please call the generic names, such as RNNOA here, and let the compiler sort out the datatypes.

[pedrohcgs]

Hi Steve,

I have just followed your advices and changed the code using the generic names. 

I don't know if this helps you cracking the problem, but if I do not set the seed nor use RNSET, the program works fine.

If I change the iseed to different values, the program freezes at different times.

Thanks 

[mecej4]

I have an observation that does not bear directly on the freeze-up problem (which I can reproduce). There is only one line where RNPOI is called. In that line, you retrieve NR (=500) random numbers, but you use only the first of those. Therefore, you could alter the call to CALL RNPOI (theta, z, NR=1).

Making this change caused the freeze-up to go away, and your program ran to a normal termination. If this is acceptable to you, you will be in a better position, although I think that the freeze-up ought to be investigated and fixed.

[pedrohcgs]

mecej4, your on-the-fly solution works very well for me!

Thanks for that!

Nonetheless, I think it is worth to investigate further the freeze-up problem.

[Steven_L_Intel1]

I am looking into that. I think it is an IMSL bug but I am not sure.

[mecej4]

The RNPOI bug seems to have been known since a few years. See

http://en.usenet.digipedia.org/thread/1128/28684/

http://coding.derkeiler.com/Archive/Fortran/comp.lang.fortran/2009-03/msg01017.html

I had completely forgotten the bug, though I had myself posted on it  in C.L.F. a few years ago.

[Steven_L_Intel1]

Well, the symptoms reported there are similar, but at the time I did not have access to the IMSL sources to diagnose it further. If I could have reproduced it I would have sent it to Visual Numerics for analysis.

Now that I do see what is going on, I can tell you that it is dependent on a uniform random number generator, called by RNPOI, returning a value extremely close to 1, in conjunction with certain values of THETA. This combination puts RNPOI into a loop. I have sent my analysis to Rogue Wave for them to investigate.

[pedrohcgs]

Steve and Mecej4, you guys rock! Thanks for all the support.

Now, let's see what Rogue Wave tell about it. Keep me posted.

Thanks again

[mecej4]

Steve,

I have found the seed and parameter values that make RNPOI freeze up before it returns a single random number. The reproducer is given below. The bug is really in IMSL code -- I can reproduce the error with a different Fortran compiler using IMSL/FNL-7.

(I also find a bug in the Forum software: I corrected a typo in the text, and that caused the line-breaks in the program source to be replaced by spaces - char(32). I had to re-paste the program source to restore the line breaks.)

[fortran] program RNPOIBUG USE RNPOI_INT USE RNSET_INT USE RNGET_INT implicit none real :: theta integer :: z(1) ! dimension(1) to match interface integer i,seed ! seed=1130536997 ! hangup sensitive to this value theta=6.84536 ! and this call RNSET(seed) do i=1,10 write(*,*)i,seed call RNPOI(theta,z) call RNGET(seed) end do end program RNPOIBUG [/fortran]

[Steven_L_Intel1]

Thanks, this is helpful! I worried that Rogue Wave might not be able to reproduce the problem if they had any small differences in the random number generator. I will pass this along.

[mecej4]

Pedrohcgs:

Intel MKL contains a Poisson RNG as part of MKL. After setting up the stream, as described in the MKL RefMan, you replace the loop (where IMSL froze up) with

[fortran] errcode=vslnewstream( stream, brng, iseed ) do j=1,nr errcode = virngpoisson( method, stream, 1, y(j), lambda(j) ) end do [/fortran]

[Steven_L_Intel1]

Rogue Wave has opened IMSL Technical Support Ticket 001-00-066751 for this. They tell me it had not been reported to them before, and that it reproduces in IMSL 7 as well.