I am interested in generating random members of two classes of $n \times n$ positive-definite matrices $A$ and $B$--the former symmetric in nature, the latter, Hermitian.

The standard (Ginibre-matrix-based) approach (Ginibre) to obtaining random realizations of $A$ is to generate an $n \times n$ matrix $a$, the $n^2$ entries of which are standard univariate normal variates. Then, the desired realization of $A$ is gotten by computing the matrix product $a a^T$, where the transpose is indicated.

To obtain random realizations of $B$, one starts with an $n \times n$ matrix $b$, the $n^2$ entries of which are of the form $c + d I$, where $c$ and $d$ are themselves standard univariate normal variates. Then, the desired realization of $B$ is obtained by computing the matrix product $b b^\dagger$, where the conjugate transpose is denoted.

(If the random realizations of $A$ and $B$ are normalized to have trace 1, then they are--in physics parlance--density matrices randomly distributed with respect to Hilbert-Schmidte measure. It is, then, of interest to test their "partial transposes" for positive-definiteness. For $n=4,6$, such tests are equivalent to ones for "separability" [PPT].)

Can such desired random realizations be accomplished with the WishartMatrixDistribution command (in conjunction with the RandomVariate command), and, if so, might there be an associated speed-up in comparison to the Ginibre-matrix-product approach outlined above?

Also, more generally, what if the classes of matrices are $k \times n$, with $k<n$? (This pertains to the PPT/separability problem for $n \times n$ matrices of rank $k$.)

  • $\begingroup$ What kind of speed gains are you looking for? The two approaches you cite for 1000 x 1000 matrices can be generated on the order of 40 ms. Did you want something significantly faster? $\endgroup$ – Carl Woll Sep 1 '20 at 18:23
  • $\begingroup$ Thanks, Carl Woll! No, the matrices in which I'm currently interested in are of much smaller dimensions--so I guess you're saying no real speed differences. Are the covariances matrices to be used in the WishartMatrixDistribution command in my application, simply unit-diagonal ones? $\endgroup$ – Paul B. Slater Sep 1 '20 at 19:51

I posed the question here to a colleague, Vahagn Abgaryan.

He responded:

"Concerning your question on WishartMatrixDistribution. First, I didn't know that such a function existed, thank you for an interesting finding. Second, this function indeed may be used to generate the real matrices. I include a small program showing the implementation of this function for full rank rebit-rebits (sorry for a haphazard and dirty code). You may see from the program that two ways are generating the same eigenvalue distributions. Interestingly enough this function seems to be 5 times slower than the direct approach (analogous the one I sent you the last time), when working with Table, and twice faster when working with ParallelTable. Anyways, I am not sure how to use it for the general complex case."

Here is the link to the illustrative Mathematica notebook that Abgaryan sent. (I put it in the Wolfram Cloud, as I was not immediately sure how to post it more directly.)



Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.