Question: Is there a parameterization of a a general qutrit 3-level system similar to: $$\rho = \begin{bmatrix} p_1 & r_{12}e^{-i\phi_{12}} & r_{13}e^{-i\phi_{13}}\\ \cdot & p_2 & r_{23}e^{-i\phi_{23}}\\ \cdot & \cdot & 1-p_1-p_2 \end{bmatrix} \tag{1.1}$$ With the constraints: $$\begin{align} & 0\leq p_j \\ & p_1+p_2\leq1 \\ & 0\leq r_{ij}\leq \sqrt{p_ip_j} \\ & 0 \leq \phi_{ij} < 2\pi \end{align} \tag{1.2}$$ ?

Today I learned a 3-level system can be parameterized, from Gell-Mann matrices, as:

$$ \rho = \begin{bmatrix} \frac13+a_3+\frac{a_8}{\sqrt 3} & a_1-ia_2 &a_4-ia_5 \\ \cdot & \frac13-a_3+\frac{a_8}{\sqrt 3} & a_6-ia_7 \\ \cdot &\cdot & \frac13-\frac{2a_8}{\sqrt 3} \end{bmatrix} \tag{2.1}$$ Subject to the constraints: $$\begin{align} & 0\leq \frac13\pm a_3+\frac{a_8}{\sqrt 3}\leq1 \\ & 0\leq \frac13-\frac{2a_8}{\sqrt 3}\leq1 \\ & \mathbf a\cdot \mathbf a\leq\frac13 \\ & 0 \leq \det(\rho) \end{align} \tag{2.2}$$

Disclaimer: I know parameterization $(1)$ is not correct, hence I am asking if there is anything similar. In particular, my problem with the Gell-Mann parameterization $(2)$ is that I have to compute a determinant each time I draw 8 random numbers... I would like something faster to generate and validate.

Background: I am doing numerical simulations, where I need to quickly generate a large number of random 3-level density operators. I was using parameterization $(2)$, with which you can easily draw the eight degrees of freedom of the qutrit.

However, the parameterization is wrong since it allows non-physical states, such as: $$\rho = \begin{bmatrix} p_1 &\sqrt{p_1 p_2} &\sqrt{p_1 p_3}\\ \cdot &p_2 &\sqrt{p_2 p_3}e^{-i\pi}\\ \cdot &\cdot &p_3 \end{bmatrix}$$ which is not positive semidefinite.

  • 2
    $\begingroup$ Do you also care about making a fair sampling? (for some definition of fair) A better starting point might be to start with a valid diagonal matrix and apply a random unitary drawn from the Haar distribution. $\endgroup$
    – DaftWullie
    Mar 15, 2023 at 12:12
  • $\begingroup$ in what sense are you saying the parameterization you gave with $p_i,r_{ij}$ is "not correct"? Looks quite correct to me. $\endgroup$
    – glS
    Mar 15, 2023 at 12:13
  • $\begingroup$ @DaftWullie thanks for the pointer. I guess I should care about fair sampling, I will look into Haar distribution. $\endgroup$
    – G Frazao
    Mar 15, 2023 at 12:18
  • $\begingroup$ @glS Thanks for the comment. I think the parameterization is not ensuring positive semidefiniteness. See the example I give in "background", it was generated using the said parameterization, but it has a negative eigenvalue. $\endgroup$
    – G Frazao
    Mar 15, 2023 at 12:21
  • 1
    $\begingroup$ ah, yes, that's right. That kind of parametrization works fine for 2x2 matrices but there's issues in higher dimensions. It's fixable using different constraints though. I actually worked that precise thing out some time ago, and if I remember correctly there is a decently nice solution. I'll try and write an answer later. Though regarding the titular question, you can take any orthonormal set of operators forming a base for the space instead of the Gellman matrices and you'll another parametrization, and the bounds are then easily writable, see physics.stackexchange.com/a/425101/58382 $\endgroup$
    – glS
    Mar 15, 2023 at 13:17

2 Answers 2


First, note that the direct answer to the question "Is there another parameterization of a qutrit 3-level system, besides Gell-Mann?" is a pretty obvious yes. There isn't anything intrinsically special about using Gell-Mann matrices to decompose qutrits, it's a basis of operators like many others (well, it has some symmetries that make it interesting in some areas, as discussed e.g. here, but this doesn't matter much for the purpose of this question). But of course, parametrising state space in dimensions larger than two is tricky due to the boundary not having a simple spherical shape, which makes the constraints on any such parametrisation a potentially delicate matter. Some of this has already been covered in What are the conditions ensuring a two-qubit density matrix is positive semidefinite?, and links therein. Here I'll try to focus a bit more on the parametrisation you mention and the associated semipositivity constraints.

Let's consider $d$-dimensional states. A very general answer to the question is that you can take an arbitrary orthonormal basis of Hermitian operators $\sigma_i\in\operatorname{Herm}(\mathbb{C}^d)$. Any state $\rho$ is then represented (or "parametrised", if you prefer) by the vector of coefficients in the decomposition $\rho=\sum_{i=1}^{d^2-1}c_i(\rho)\sigma_i$. Note that here the "orthonormality" means we're requiring $\operatorname{tr}(\sigma_i\sigma_j)=\delta_{ij}$. It is also common to use operators normalised as $\operatorname{tr}(\sigma_i\sigma_j)=d\delta_{ij}$, but this doesn't change much in the discussion, so I'll stick with actual orthonormal operator bases here.

Given the normalisation constraint on density matrices, it is always convenient to include the identity operator in the basis. Rescaling to ensure normalisation, we can thus assume without (much) loss of generality that $\sigma_0=I/\sqrt d$, and then orthogonality implies that all $\sigma_i$ with $i>0$ are traceless. Then, for any such basis, you can characterise the boundary of state space in "hyperspherical coordinates" with the constraint $$\|\vec c(\rho)\| \le \frac{1}{|\min_k \lambda_k(\rho)|},$$ where $\lambda_k(\rho)$ is the $k$-th eigenvalue of $\rho$. You can find a proof of this e.g. here.

Regarding a parametrisation like your (1.1), in terms of probabilities and phases, consider rewriting it as $$\rho = \begin{bmatrix} p_1 & \sqrt{p_1 p_2} r_{12}e^{-i\phi_{12}} & \sqrt{p_1 p_3 }r_{13}e^{-i\phi_{13}}\\ \cdot & p_2 & \sqrt{p_2 p_3}r_{23}e^{-i\phi_{23}}\\ \cdot & \cdot & 1-p_1-p_2 \end{bmatrix}, $$ where all parameters are real. Sylvester's criterion now tells you that for positive semidefiniteness is sufficient to check all the principal minors. The conditions associated to 1x1 submatrices amount to $(p_1,p_2,1-p_1-p_2)$ being a probability vector, that is, $0\le p_1,p_2,p_1+p_2\le 1$. The conditions associated to 2x2 submatrices give you $0\le r_{ij}\le 1$. Finally, the condition on the full determinant can be put in the form $$\|\vec r\|^2 \le 1 + 2r_{12} r_{13} r_{23} \cos(\phi_{13}-\phi_{12}-\phi_{23}).$$ Here $\|\vec r\|^2\equiv r_{12}^2+r_{23}^2+r_{13}^2$. Note how this choice of parameterization with the $\sqrt{p_i p_j}$ factors gives conditions on the rest of the parameters which are independent on the probabilities themselves. This kind of thing probably works also in higher dimensions, and maybe you can write the associated general expressions, but it might get quite a bit trickier, and I'd guess you'll have to keep track of how different permutations classes contribute to the different determinants. I haven't actually tried.

Another advantage of this type of representation is that it gives a particularly simple expression for the purity: $$\operatorname{tr}(\rho^2) = \sum_{j=1}^3 p_j^2 + 2\sum_{j<k} p_j p_k r_{jk}^2 = 1 + 2\sum_{j<k} p_j p_k (r_{jk}^2-1),$$ and thus the phases $\phi_{jk}$ do not affect the purity at all, while $r_{jk}=1$ ensures pure states, because you can then complete the square.

  • $\begingroup$ Thanks a lot. That is precisely the kind of parameterization I was looking for. And $\|\vec r\|^2 \le 1 + 2r_{12} r_{13} r_{23} \cos(\phi_{13}-\phi_{12}-\phi_{23})$ sure is a better-looking constraint, and with more insight than $\text{det}(\rho)>0$. Note however that the phases $\phi_{jk}$ do affect the purity! The necessary and sufficient condition for purity is $r_{jk}=1$ and $(\phi_{13}-\phi_{12}-\phi_{23})= 0$, otherwise the parameterization fails the determinant constraint. $\endgroup$
    – G Frazao
    Mar 16, 2023 at 16:04
  • $\begingroup$ $(\phi_{13}-\phi_{12}-\phi_{23})= 0 \mod 2\pi$ $\endgroup$
    – G Frazao
    Mar 16, 2023 at 16:11
  • 1
    $\begingroup$ @GFrazao eh, not quite. The phases $\phi_{ij}$ never affect the "purity", that is, the quantity $\operatorname{tr}(\rho^2)$. But they do certainly affect whether $\rho$ is positive semidefinite. As you can easily verify numerically, changing values of $\phi_{ij}$ outside that condition you get non-positive-semidefinite matrices, which however still verify $\operatorname{tr}(\rho)=\operatorname{tr}(\rho^2)=1$. $\endgroup$
    – glS
    Mar 16, 2023 at 19:59

Updated answer after DaftWullie's comment:

I think I over-complicated; if $U$ is a random unitary in $SU(3)$, $U\vert0\rangle$ is good enough to get a random pure state, so you could get a random pure density matrix with

$$\rho = U\vert0\rangle\langle0\vert U^\dagger$$

and then convexly combine pure density matrices to get mixed ones.

On the other hand, what DaftWullie was really suggesting seems even less complicated, and so more reasonable; just pick a random basis that diagonalizes the matrix.

$$\rho = U (a \vert{0}\rangle\langle{0}\vert + b\vert{1}\rangle\langle{1}\vert + (1-a-b)\vert{2}\rangle\langle{2}\vert) U^\dagger$$

and the comments below about the parametrization of random unitaries still apply.

Old answer:

I don't have enough reputation to comment, unfortunately, so consider what's below as an answer stub, rather than a fully fledged answer, though I think it might be a good first step.

As DaftWullie mentioned, you can get a pure density matrix more-or-less easily, by choosing suitable amplitudes in a random basis; thus if you have the ability to generate a Haar-random unitary matrix $U$, you can do

$$ \rho = U (a \vert0\rangle + be^{i\phi_1}\vert1\rangle + e^{i\phi_2}\sqrt{1-a^2-b^2}\vert2\rangle) \cdot \mbox{c.c.} $$

for $a,b \in [-1,1]$ satisfying $a^2 + b^2 \leq 1$.

Parametrizing the space of unitaries seems doable, but not easy; see this answer.

Once you can generate random pure density matrices, you can make use of the fact that the density matrix space is convex; choose two matrices and interpolate between them. I doubt this will give you a uniform distribution though, simply based on the intuition that getting the fully mixed matrix will require a very specific choice of vertices for the interpolation.

I hope someone can build on this and help you further.

  • $\begingroup$ Did you really mean this form of $\rho$? I was thinking $U(a|0\rangle\langle 0|+b|1\rangle\langle 1|+(1-a-b)|2\rangle\langle 2|)U^\dagger$ for $0\leq a,b$ and $a+b\leq 1$. $\endgroup$
    – DaftWullie
    Mar 15, 2023 at 13:38
  • $\begingroup$ I did mean this form of $\rho$, yes (well, sort of, I'd forgot the phases and had a wrong range for $a$,$b$, which I've fixed now), since I was thinking about pure states, but I see now what you mean now (form a direct diagonalization point of view). $\endgroup$
    – MikeEVMM
    Mar 15, 2023 at 14:15
  • $\begingroup$ Thanks! However... 1) the linked answer points to a report on how to uniformly generate a random unitary, not to a parameterization. 2) I really need a parameterization so I can explore the parameter space with an optimization algorithm. $\endgroup$
    – G Frazao
    Mar 15, 2023 at 16:35
  • $\begingroup$ Are you sure? I'm asking because of the paragraph at the end of section 2.3 of Maris Ozols's document, linked in the first answer to that question. Since you can parametrize 2×2 unitaries, maybe ref. 8 will be of help? ar-tiste.com/cybenko2001.pdf Likewise ref. 6 can be found at chaos.if.uj.edu.pl/~karol/pdf/ZK94.pdf $\endgroup$
    – MikeEVMM
    Mar 15, 2023 at 17:26
  • $\begingroup$ I am not sure :-) I also found this paper which does exactly what @DaftWullie and you propose. I will read it tomorrow and update later. Thanks for all the effort. $\endgroup$
    – G Frazao
    Mar 15, 2023 at 19:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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