Re: Two suggestions...

[James Theiler <jt at lanl dot gov>]

On Tue, 10 Jun 2003, Robert G. Brown wrote:
] A couple of consistency-check questions or suggestions as I'm writing a
] pretty complicated GSL-based numerical program.
]
]<...>
]
]   b) It is numerically much more efficient to evaluate random numbers a
] vector at a time.  This saves two ways.  By creating an "empty" random
] number generator in the gsl format and adding it to my test/benchmark
] program I measure overhead on the function call alone to be between 10
] and 30 nanoseconds on CPUs that clock from 800 MHz to >2 GHz, and this
] can is anywhere from 1/3 to 2/3 of the total time required to evaluate
] e.g. mt19937.  Simply using one function call to generate a full page of
] rands at a time would thus result in a speedup by as much as 3.  Next,
] many rng's can use pipelines and cache much better in a vector loop
] where the CPU, memory type, and compiler permit, especially if they
] don't have to give up CPU registers in between.  This could speed up
] rng's by another factor of 2, for as much as 6 overall -- from perhaps
] 30 nanoseconds (or more) to as little as 5 or 10 nanoseconds.  This kind
] of speedup matters in lots of Monte Carlo applications.
]
] Is there any thought about creating a vector-rng interface with perhaps
] an inline conditional macro instead of a subroutine call?
]


On 11 Jun 2003, Gerard Jungman wrote:
] <...>
] Also, I thought the comment about vectors of random numbers
] was probably correct; it might be a bit tedious to add a
] vector version of the rng functions, but it seems like
] that's the least we can do, given all the good work that
] has already gone into that library component.

Hi,
  For what it's worth, I'm inclined to agree with Gerard here, though
I am sympathetic to the argument that we should be reluctant to add
new functions to what is already a pretty extensive API. The case for
extending the API could be made stronger by demonstrating an actual
speedup.  I think Robert's argument above is credible, but at this
point somewhat theoretical?  I might add that batch computation of
random numbers in the randist library should also be considered.  I
know from the binomial distribution I did, and the Poisson that I have
been promising to do[*], that substantial speedups can be achieved
this way.
  By the way, that last factor of two could still be achieved using
existing GSL, and looping over the elements of a vector to create the
random numbers, and then letting vectorization do its stuff when it
comes time to use those random numbers. But I have to agree that it
would be more elegant (and up to 3x faster, it sounds like) to produce
the vector all at once.
  Another question I guess is how hard it would be to make the change
(partly it's a question of how tedious it would be, but mostly the
issue would be how maintainable the code would be afterward).  Perhaps
a trial implementation would address this issue as well.

regards,
jt

[*] essentially done, but want to do some cleanup before I contribute
it, and -- surprise -- other things have come up in the last few weeks
that have kept me from finally finishing it.

---------------------------------------------
James Theiler                     jt@lanl.gov
MS-B244, NIS-2, LANL        tel: 505/665-5682
Los Alamos, NM 87545        fax: 505/665-4414
----- Space and Remote Sensing Sciences -----