How to construct symmetric and positive definite $A,B,C$ such that $A+B+C=I$? In an attempt to formulate a answer to this (in)famous question

How does one prove the matrix inequality $\det\left(6(A^3+B^3+C^3)+I_{n}\right)\ge 5^n\det(A^2+B^2+C^2)$?

I'm trying to construct three $n\times n$ matrices $A,B,C$ that are (a) symmetric, (b) positive definite, (c) add to $I_n$ .
Note that I've already decided to restrict attention to the reals and I have replaced Hermitian by symmetric (which IMO is
difficult enough).

My unsuccessful tries are a wild mixture of two extremes:
Make a random square matrix $H$ and form $A = H^TH$ . Make another random square matrix $H$ and form $B = H^TH$ .
In the same way, form $C = H^TH$ . Then $A,B,C$ are symmetric and positive definite. But in general $A+B+C \ne I$ .
Generate random numbers for $A_{ij} = A_{ji}$ , $B_{ij} = B_{ji}$ and form $C_{ij} = C_{ji} = I_{ij}-A_{ij}-B_{ij}$ .
Then $A,B,C$ are symmetric and $A+B+C = I$ , but it cannot be guaranteed that these are positive definite matrices.
So the question is: how can the three requirements (a) , (b) , (c) be fulfilled at the same time, while keeping $A,B,C$ yet as
random as possible? My plan is to do numerical experiments and eventually find a counter example. I have all the ingredients to do it, except this.
 A: This is easy. Denote by $\operatorname{diag}(x)$ the diagonal matrix whose diagonal is the vector $x$. Assuming that $A$ is diagonal, every valid triple $(A,B,C)$ can be expressed as follows:
\begin{align*}
A&=\operatorname{diag}(a),\\
R&=(I-A)^{1/2},\\
B&=RU\operatorname{diag}(b)U^\ast R,\\
C&=I-A-B,
\end{align*}
where $a$ and $b$ are two vectors whose entries lie inside $(0,1)$ and $U$ is a unitary matrix. Note that $C=RU\left(I-\operatorname{diag}(b)\right)U^\ast R$. Hence it is positive definite.
So, to generate a random sample, all you need are random instances of $a,b$ and $U$. The unitary matrix $U$ can be obtained by performing a singular value decomposition on a random square matrix, or you may generate it using the method mentioned in Wikipedia.
A: Another option the following.

*

*Pick arbitrary random matrices $X$, $Y$, $Z$ (either real or complex depending on whether symmetric or Hermitian matrices are needed).

*Set $A_0 = X^\dagger X$, $B_0 = Y^\dagger Y$, $C_0 = Z^\dagger Z$. These are
positive Hermitian, but do not necessarily add up to one.

*Compute $D = A_0 + B_0 + C_0$ and set $A = D^{-1/2} A_0 D^{-1/2}$, $B = D^{-1/2} B_0 D^{-1/2}$, $C = D^{-1/2} C_0 D^{-1/2}$.

By construction $A+B+C = D^{-1/2} D D^{-1/2} = I$.
The minor downside of this construction (depending on the distribution used to pick $X$, $Y$, $Z$) is that it is theoretically possible that $D$ has a 0 eigenvalue. That seems unlikely for most of the distributions for $X$, $Y$, $Z$, since that would imply that the corresponding eigenvector is
the eigenvector with 0 eigenvalue for all 3 of $X$, $Y$, $Z$ simultaneously. If that happens, one may restart the sampling from the beginning, or add a small positive number to $A_0$, $B_0$, or $C_0$.
