Algorithm for generating positive semidefinite matrices I'm looking for an efficient algorithm to generate large positive semidefinite matrices.
One possible way I know of is:


*

*generate a random square matrix

*multiply it with its transpose.

*calculate all eigenvalues of the result matrix and check if all of them are non-negative.


However, this approach is infeasible given a large matrix, say $1000 \times 1000$ or more. 
Could anyone please suggest an efficient way to generate a positive semidefinite matrix?
Is there any MATLAB function for this job?
Thanks,
 A: *

*Generate a random symmetric matrix, determine eigenvalue bounds via, say, Gerschgorin, and then shift the diagonal elements by an appropriate amount determined from the Gerschgorin bound of the leftmost eigenvalue.

*Generate a diagonal matrix with random nonnegative elements from a distribution of your choice, and perform a similarity transformation with a Haar-distributed pseudorandom orthogonal matrix.

*Generate a diagonal matrix with random nonnegative elements from a distribution of your choice, and perform a few sweeps of the (cyclic) Jacobi algorithm, with randomly generated rotation matrices $\begin{pmatrix}c&-s\\s&c\end{pmatrix}$ (e.g., randomly generate a $c\in [-1,1]$ and calculate a corresponding $s$ through $c^2+s^2=1$).
A: Generate some "random" vectors $\mathbf v_1,\dots, \mathbf v_m$ and "random" non-negative scalars $c_1, \dots, c_m$ and compute 
$$\mathbf P=c_1 \mathbf v_1\mathbf v_1^\top+\cdots+c_m \mathbf v_m\mathbf v_m^\top$$
A: Pick an inner product in $\mathbb R^m$ or in $\mathbb C^m$, a set of vectors $v_1$, $\dots$, $v_n$ in that space, and consider the $n\times n$ matrix $A=(a_{i,j})$ with $a_{i,j}=\langle v_i,v_j\rangle$. This is called the Gramian matrix of the vectors you started with, it is always positive semidefinite, and in fact every positive semidefinite matrix is the Gramian matrix of some set of vectors.
A: Any matrix multiplied by it's transpose is going to be PSD; you don't have to check it.  On my computer raw Octave, without SSE, takes 2 seconds to multiply a 1000x1000 matrix with itself.  So not all that infeasible.
If you don't like that, you can always just generate a random diagonal matrix.  That's sort of the trivial way, though :)  What do you need the matrix for?  If it's as test input to another algorithm, I'd just spend some time generating random PSD matrices using the above matrix-matrix multiplication and save the results off to disk.
