How to approx. decompose a sym. p.d. matrix M into X'X? M: pxp  symmetric p.d. matrix with unit diagonals
n:  number much smaller than p
Want a nonrandom nxp matrix X such that X'X is 
close to M element-wise.  If n gets larger, hopefully 
difference gets smaller.  
I have a no so good method. Get the eigendecomposition 
M = D'VD, and rearrange the diagonals of V in decreasing order 
and so the rows of D. Take only the first n  rows of D
as X. This method works only if the eigenvalues decrease 
very fast.
Any input is appreciated. Thanks a lot! 
 A: Unfortunately, you can't do any better than your "no so good method."  It's a standard result that this is the best (as measured by the Frobenius norm) rank $n$ approximation to $M$. 
A: Not a new solution but an elaboration from a computational viewpoint.
From a practical point of view one can solve your problem as follows:
Firstly, we do not care about the diagonal values being $1$.
Secondly, we enlarge the problem to real symmetric matrices with 
no strictly negative eigenvalues.
Such a matrix $A$ of size $p\times p$ is then given by $A=\sum_i \lambda_i v_i^tv_i$
where $v_1,\dots,v_p$ is an orthogonal basis of eigenvectors with eigenvalues $\lambda_1\geq \lambda_2\geq\dots$
of $A$.
Suppose there is a fast way for computing an eigenvector $v_1$ (of norm $1$) 
associated to the largest eigenvalue $\lambda_1$ of $A$. (One can for instance consider the projective
limit of $A^k v$ for $v\in\mathbb R^p$ a generic vector.)
The largest eigenvalue $\lambda_1$ is then given by 
$$\lambda_1=\sum_{i,j}A_{i,j}{P(1)}_{i,j}=\sum_{i,j}A_{i,j}{(v_1)}_i{(v_1)}_j$$
where $P(1)=v_1^t v_1$ is the orthogonal rank one projector associated to $v_1$.
Replacing $A$ by $A-\lambda_1P_1$ and iterating one gets $P_2,P_3,\dots$ and $\lambda_2,\lambda_3,\dots$.
Your solution is then given by $\sum_{i=1}^p\lambda_iP(i)=\sum_{i=1}^p\lambda_iv_i^tv_i$.
