@acer : Indeed, all of what I discussed is specific to user-defined distributions, but your question applies to both the predefined and the custom distributions, in principle: for both, there are situations where there is data we want to keep between invocations of the C level code. However, we currently have facilities for doing so only for the predefined distributions. Fixing this has been on my wish list for a while, but it hasn't bubbled to the top quite yet.
So in the mean time, for custom distributions, it is still best to, as you say, create a sample that is as big as your memory can bear and then slice and dice that to obtain your subsamples. For predefined distributions, it should be almost as efficient if you use this underdocumented calling sequence:
X := Statistics[RandomVariable](SomePredefinedDistribution(alpha, beta));
p := Statistics[Sample](X); # This does most of the library-level precomputation
for i to 10^5 do
sample := p(10^4):
This took 65s for SomePredefinedDistribution(alpha, beta) = Normal(0, 1) on one machine where I tested it (with the process statement omitted), versus 63s for 100 samples of size 10^7. We will try to document this calling sequence (where Sample returns a procedure) better in a future version of Maple.
Hope this helps,