Calling random_number ( ) in pure procedure

I need call random_number (r) in my pure subroutine. r is an array of 2x2. But I get this error.

Error: Subroutine call to intrinsic 'random_number' at (1) is not PURE

One way I find is to pass this array into the pure procedure and declare it as intent in. That works fine. But is there any way to call it in the pure procedure?

You cannot call the random_number intrinsic in a pure procedure, but you could write a pure RNG that is slightly modified from a non-pure one, as discussed in Pure random number generators, and call that.

Alternatively, since a pure procedure can access a module variable but not change it, you could create a wrapper that calls random_number, sets a module variable, and then calls the pure subroutine, as shown below.

module m
implicit none
real :: xran
contains
subroutine wrap_add(x)
real, intent(in out) :: x
call random_number(xran)
print*,"xran =",xran
call add(x)
end subroutine wrap_add
!
pure subroutine add(x)
real, intent(in out) :: x
x = x + xran
end subroutine add
end module m
!
program main
use m
implicit none
real :: x
x = 10.0
call wrap_add(x)
print*,"x =",x
call wrap_add(x)
print*,"x =",x
end program main

sample output:

 xran =  0.349059999    
 x =   10.3490601    
 xran =  0.555532038    
 x =   10.9045925
1 Like

This solution is not thread-safe, which may be needed, depending on OP’s goals. There is !$omp threadprivate, but OpenMP trainers recommend against using it if avoidable.

1 Like

If you really want a pure procedure, then that would be the standard way to work around the limitation.

The random_number() intrinsic saves an internal state, so it cannot be pure, and thus it cannot be called from a pure procedure.

1 Like

If all you want is the procedure to be elemental just make it impure and it should work fine with random nuumber in it. Otherwise I would use a pure RNG, as @Beliavsky recommended.

1 Like

This is true for all random number generators of which I’m aware - it’s not restricted to the intrinsic procedure

If the PURE procedure requires many random deviates another approach is to initialize the RNG each time the procedure is called, either from data or (better) from a custom-built seed - this essentially initializes the RNG state and throws it away during each function call. This is inefficient if only a few random numbers are needed.

As an alternative, you could call a custom PRNG that expects the seed to be passed as a dummy argument with each call. You would then need to provide a seed array to your pure elemental routine.

1 Like

Yes, this would be true for any pseudorandom number (PRNG) procedure.

However, I’m curious about a random number generator (RNG) maintained by the processor based on, for example, thermal noise, or some kind of nuclear decay mechanism, or cosmic ray detection. That kind of RNG does not really have an internal state that is saved. Its bits are continually being updated, and (I’m not sure about this) copying those bits at any moment in time does not affect their later values. Would it be legal for a pure fortran procedure to access the value of such a RNG?

1 Like

That is exactly what I have been doing.

We need to introduce thermal fluctuations to trigger nucleation events. This is very crucial, especially in the case of secondary and tertiary dendritic arm formation in Phase-field simulation.

https://www.sciencedirect.com/science/article/abs/pii/S1003632610602087?via%3Dihub

But here I am more interested in the code structure and readability. I find it pleasant if the implementation details and extra procedure calls are hidden in the main program. Also, it avoids passing on the additional argument. But here it seems inevitable!

If it is any useful, if you have enough memory, you can try to generate a sufficiently big number of prn, store them in an array, or a file, and take them from this storage. This way you can reject the correlated numbers, store only the uncorrelated ones.

Depending on your particular circunstances this may viable or useful, or maybe not. In any case, taking the numbers from an array is faster than calling a function.

1 Like

The array suggestion will work. The original poster already did that with a 2x2 array argument. However, I think the other option of writing the values to a file will not work. I think that reading from an external file counts as a side effect that is not allowed within a pure procedure. The side effect is not reading the values but rather positioning the file.

My previous question involved a RNG that is maintained by the processor but does involve saving a state. Here is an article from intel about this. Behind Intel’s New Random-Number Generator - IEEE Spectrum Other types of hardware has been used for this purpose in the past too. This basically defines a global register that can be read, but not modified, where the bits are random (or as random as practically possible). A pure procedure is allowed to read (but not modify) global entities, such as module variables defined within the fortran standard. My question is whether a pure procedure in fortran would be allowed to read this thermal noise register? Of course, the compiler would need to package it somehow, for example as a VOLATILE variable within an intrinsic module, to make it accessible through fortran.