multi-variate distribution

[Ulrich Drepper <drepper at gmail dot com>]

The distribution functions included so far are only the basic ones.
For simulations what is needed are multiple variables and they usually
are correlated.  This is what multi-variate distributions are about.
The most commonly used is the multi-variate normal distribution.  The
attached patch implements such a distribution.

The details for the implementation differs from that of other
distributions.  The main problem is that the produced value are
(except in the degenerated case) not a single value.  Given that the
operator[] has to be there and it has to return all the results the
result_type type has to be large enough.  I chose std::array for that.
 This means min() and max() return vectors, too, even though all
dimensions use the same type and have the same limits.

Then there are the constructors.  There are two different parameters:
a k-dimensional mean vector and the description of a kxk covariance
matrix.  There are some other constructors in <random> for more
complex parameters but then there is only one.  What to do for two
complex parameters?  In the patch so far I have used the same method
for the vector and the matrix.  Should there be a mix as well?

It is usual that the covariance matrix is specified in three different
ways: a full matrix, an upper/lower triangle matrix, or just the
diagonal.  In the implementation the type of the initializer is
deduced from the number of values passed.  This always works, even for
the degenerated case k=1.  Is this something which fits in with the
general design of the library?  The alternative would be to have a
separate parameter which is just an enum.  I thought this to be
inelegant.

Another problem is the the constructor can fail.  The matrix must be
positive-definite.  There are other functions in <random> which fail
and just throw a runtime_exception.  Is this the right way?

The normal_mv_distribution type has two template parameters.  The
first is the dimension.  The order differs from that of std:array.
Not so nice but the dimension cannot have a default parameter.  It
would be possible to have just one parameter.  It's possible to
rewrite the code to have the dimension be part of the param_type.
This would make some code slower.  The computations in operator[] and
__generate are likely to be fully optimized and linearized for small
dimensions.  What is the recommendation here?

Dropping the dimension makes it easier to write optimized __generate
functions because the base type of the iterator would be _RealType* or
something like this.  But it means the type differs from result_type
which in this case would have to be std::vector or something like
this.  This is costly and which is mainly why I chose to use a
template parameter for the dimension.  Performance is chief
consideration as far as I am concerned.

One could argue about the parameter to use.  The covariance matrix is
probably the most common.  The only argue to not use it that I can
imagine is that it limits the parameters to those that are valid
(positive-definite matrix).  I cannot think about a valid reason why
someone would want to use a parameter that cannot be expressed with a
covariance matrix.  The nice separation between creating parameters
and the performing calculations make it a no-problem that the
covariance matrix has to be transformed first.  It's dealt with
efficiently.

Any opinions?  Otherwise I'll submit the patch for inclusion.

Attachment: d-random-normmv
Description: Binary data