Re: cellfun and parcellfun

[Robert I A Patterson]

David,

On Thu, 2010-06-10 at 21:52 +0200, David Bateman wrote:
> Jaroslav Hajek wrote:
> > On Thu, Jun 10, 2010 at 7:24 PM, Robert Patterson <address@hidden> wrote:
> >   
> >> Hi,
> >>
> >> I am writing some Monte Carlo simulation code.  The problem can be
> >> broken down via domain decomposition (several different sub-problems
> >> that only need to communicate occasionally).  I am implementing this
> >> decomposition via cellfun, in order to avoid writing explicit loops.  I
> >> hope to parallelise with parcellfun from the general package on
> >> octaveforge.
> >>
> >> Assuming I can reseed the random number generator with the same values
> >> each time I run a calculation:
> >>
> >> 1. Will I get the same random number at each point in my code when using
> >> cellfun as I would with a simple loop, or is the order at least stable
> >> on any one implementation?
> >>
> >>     
> >
> > Yes. Cellfun evaluates sequentially upwards, just like a loop.
> >
> >   
> >> 2. Is it safe to use parcellfun with the built in random number
> >> generator?
> >>
> >>     
> >
> > Yes. parcellfun is not based on threads, but on process cloning
> > (fork). Using any function is therefore safe. However, all cloned
> > process will experience the same sequence of random numbers. The
> > easiest way to avoid this is to generate all number beforehand or
> > generate a seed for each input. Eventually, I could create a more
> > general approach to solve this.
> >
> >   
> If you call parcellfun once and before a call to rand, the state of the 
> random generator will be undefined.. The first call to rand in the 
> forked process should end up with a different state vector, either based 
> on /dev/urandom if that entropy pool is available, or the lsb of the 
> clock, which should still have enough difference too ensure that the 
> forked processes have different random sequences.
> 
> It might be nice to add a feature like rand("state", "reset") or 
> rand("state", []) that would reinitialize the state of the generator so 
> that we could ensure that state is different in each forked process even 
> if rand had been called previously
> 
> D.
> 
I understand the reasons, although it leads to quite awkward behaviour
from a Monte Carlo point of view. Imagine I run 8 sample paths using 8
processes via parcellfunc.  If I have already used some random numbers,
all 8 paths will be identical, but if this happens to be the first thing
I do when I start Octave I will get 8 different samples (provided each
process gets different data from /dev/urandom).  I wonder if some kind
of note in the parcellfun documentation might be appropriate:
"Note that each process will inherit a copy of the random number
generator state that existed before the call to parcellfun.  If that
state is uninitialised, each process will perform the default
initialisation of its random number generator separately.  In the latter
case it is possible, but not certain that each process will initialise
its random number generator differently."

Robert