Improper density matrix with negative eigenvalues

[Haorong Wu]

TL;DR Summary

I have a matrix with negative eigenvalues. This makes it inproper for a density matrix. Hence I am trying to correct it.

Hi, there. I am working with a model, in which the dimension of the Hilbert space is infinite. But Since only several states are directly coupled to the initial state and the coupling strength are weak, then I only consider a subspace spanned by these states.

The calculation shows that the trace of the matrix remains 1 which is good for a density matrix. However, the final matrix has a negative eigenvalue, -0.0661. Other eigenvalues are 0.5496, 0.4671, 0.0493. This means that the matrix will not be semi-positive. Hence it can not be a density matrix.

I guess this is due to the fact that I only consider a small fraction of the full-space.

Then in order to construct a density matrix, I try the following procedure. First, I eigendecompose the matrix. Second, I throw away the eigenvector with negative eigenvalue since the negative eigenvalue is small. Then, I reconstruct a matrix by using the remaining eigenvectors and eigenvalues. Finally, I normalize the matrix by imposing the condition $\rm{tr} (\rho)=1$.

This makes it to be a density matrix. However, I am not sure whether this procedure reasonable.

 

[stevendaryl]

Summary:: I have a matrix with negative eigenvalues. This makes it inproper for a density matrix. Hence I am trying to correct it.

Hi, there. I am working with a model, in which the dimension of the Hilbert space is infinite. But Since only several states are directly coupled to the initial state and the coupling strength are weak, then I only consider a subspace spanned by these states.

The calculation shows that the trace of the matrix remains 1 which is good for a density matrix. However, the final matrix has a negative eigenvalue, -0.0661. Other eigenvalues are 0.5496, 0.4671, 0.0493. This means that the matrix will not be semi-positive. Hence it can not be a density matrix.

I guess this is due to the fact that I only consider a small fraction of the full-space.

Then in order to construct a density matrix, I try the following procedure. First, I eigendecompose the matrix. Second, I throw away the eigenvector with negative eigenvalue since the negative eigenvalue is small. Then, I reconstruct a matrix by using the remaining eigenvectors and eigenvalues. Finally, I normalize the matrix by imposing the condition $\rm{tr} (\rho)=1$.

This makes it to be a density matrix. However, I am not sure whether this procedure reasonable.

I would want to know where the negative eigenvalue came from. Is that some error of approximation, or what?

 

[Haorong Wu]

I would want to know where the negative eigenvalue came from. Is that some error of approximation, or what?

I would like to know that too. I am working on equations of motion for density matrix. I have work them out. These are several coupled ODEs. Then by numerical solving the ODE, MATLAB gives me the final matrix with a negative eigenvalue.

Possible sources for the negative eigenvalue are:
a) in my model, there are some approximation. Some higher terms of perturbation have been ignored;
b) instead of the infinite Hilber space, I only consider a six dimensional subspace;
c) maybe the computational error caused by the MATLAB algorithm, i.e., ode45.

 

[PeterDonis]

I have a matrix with negative eigenvalues. This makes it inproper for a density matrix. Hence I am trying to correct it.

I'm not sure how we can help if you don't show us the matrix and how you got it.

 

[vanhees71]

I would like to know that too. I am working on equations of motion for density matrix. I have work them out. These are several coupled ODEs. Then by numerical solving the ODE, MATLAB gives me the final matrix with a negative eigenvalue.

Possible sources for the negative eigenvalue are:
a) in my model, there are some approximation. Some higher terms of perturbation have been ignored;
b) instead of the infinite Hilber space, I only consider a six dimensional subspace;
c) maybe the computational error caused by the MATLAB algorithm, i.e., ode45.

Is your time evolution unitary? If you simply cut the time evolution to a tiny subspace it's pretty probable that you violate unitarity, and then of course you are likely getting a wrong spectrum for the time-evolved stat. op.

Of course, we cannot really help without having a proper description of the problem you are dealing with.

 

[Haorong Wu]

Hi,@stevendaryl , @PeterDonis , @vanhees71.

I could not edit the post now. So I just rewrite it here.

I am studying the propagation of light in some medium. The state of light can be $\left | m \right >$ where $m$ stands for some OAM basis ${l_m,p_m}$.

My final result of the euqations of motion is given by $$\partial_z \rho_{mn}= L_{ns}\rho_{ms}+L_{ms}^*\rho_{sn} $$ where Einstein summation convention is used, and $L_{ns}$is some complex function depending on $z$.

Should I write how I derive the equation? But that would be very long.

I find that $L_{ns}$ could be nonzero only if $|n-s|\le 1$. Therefore, although OAM basis is infinite, I can only consider $L_{n,n-1}$,$L_{n,n}$,$L_{n,n+1}$, etc.

If the initial state is in $\left |0 \right >$. (OAM basis need two indices, but I will write only one of them for simplicity) Then I could only consider the space spanned by $\left |0 \right >$,$\left |1 \right >$. But I am afraid I may ignore too much part, so I also include $\left |2 \right >$

Then the equations of motion will be nine coupled ODE.

I solved them in MATLAB with ode45.

The initial density matrix is given by
$$\rho(0)= \begin{pmatrix}
1&0&0 \\
0&0&0 \\
0&0&0
\end{pmatrix} .$$

The final density matrix is calculated as
$$\rho(z_0)= \begin{pmatrix}
0.767 &0.118 - 0.186i &-0.123 - 0.238i\\
0.118+ 0.186i& 0.169 & 0.031 - 0.098i\\
-0.123 + 0.238i &0.031 + 0.098i &0.037
\end{pmatrix} .$$

$\rm{tr}\rho(z_0)<1$ means some information leaks to other states.

The eigenvalues of the final matrix are 0.9314, 0.1004, and -0.0583. The last one cause the matrix to be inproer for a density matrix.

Also, $L_{ns}$ is not unitary. It can be very small, or quite large, depending on parameters $z$, $w_0$, $k$, etc.

 

[A. Neumaier]

Solve instead the 3x3 system of ODEs
$$\partial_z \psi_{mn}=L_{ms}^*\psi_{sn}$$
Then $$\rho=\psi\psi^*$$ satisfies your equation and is Hermitian positive semidefinite.

 

[Haorong Wu]

Thanks, @A. Neumaier .

This solve the semidefinite problem. However, the main aim of my study is the decoherence of the system. This can not be revealed by the evolution of states. Hence I have to stick with the evolution of density matrix. Sad.

 

[Twigg]

If I'm not mistaken, the equations of motion $$\partial_z \rho_{mn} = L_{ns} \rho_{ms} + L_{ms}^*\rho_{sn}$$ do not allow unitary evolution (makes sense, given that you're interested in decoherence).

Do a quick check and see what the eigenvector corresponding to the negative eigenvalue is. I'll bet it's something of the form (small,smaller,big), just eyeballing your final matrix. I would interpret this to mean that your procedure might be (I ain't betting on it) a good approximation for $\rho_{00}$, $\rho_{01}$, and $\rho_{11}$, but it's a terrible approximation on any matrix element with a 2 in it.

If you want to check this, try incrementally raising the dimension on the matrix and seeing if $\rho_{00}$ converges. Hopefully, the negative eigenvalue will converge to a non-negative eigenvalue. If not, it might be worth investigating L more carefully.

 

[A. Neumaier]

Thanks, @A. Neumaier .

This solve the semidefinite problem. However, the main aim of my study is the decoherence of the system. This can not be revealed by the evolution of states. Hence I have to stick with the evolution of density matrix. Sad.

??

You get the evolution of the density matrix at any time by transforming the evolution of $\psi$ at this time.

 

[Twigg]

@A. Neumaier I believe the OP is saying that the system is expected to evolve into a mixed state and $\psi$ is not defined. Am I misunderstanding something here?

 

[Haorong Wu]

@A. Neumaier When I first studied decoherence, I tried to evolve a state and then calculate the density matrix. However, when a state $\left | \psi \right >$ evolves into another state$\left | \psi ' \right >$, the density matrix $\left | \psi ' \right >\left < \psi ' \right |$ will always be a pure state, as @Twigg pointed out. No decoherence can be deduced from them.

If I want to calculate decoherence by using the evolution of states, then the state must be a combination of system and environment, and at last, the environment part must be traced out.

Unfortunately, it is difficult to write a state for environment in my case.

@Twigg Thanks! I think you are right. The two positive eigenvalues correspond to two states with highest probability amplitude, about 0.9, in state $\left | 0 \right >$ and $\left | 1 \right >$, while the negative one corresponds to state $\left |2 \right >$. I will try increase the dimension.

 

[Twigg]

Sounds good! Hope it goes well!

By the way, in your equation of motion, does "*" mean "complex conjugate" or "conjugate transpose"?

 

[Haorong Wu]

Sounds good! Hope it goes well!

By the way, in your equation of motion, does "*" mean "complex conjugate" or "conjugate transpose"?

Thanks! It means complex conjugate. I would use "$\dagger$" to indicate conjugate transpose.

 

[vanhees71]

If I'm not mistaken, the equations of motion $$\partial_z \rho_{mn} = L_{ns} \rho_{ms} + L_{ms}^*\rho_{sn}$$ do not allow unitary evolution (makes sense, given that you're interested in decoherence).

Do a quick check and see what the eigenvector corresponding to the negative eigenvalue is. I'll bet it's something of the form (small,smaller,big), just eyeballing your final matrix. I would interpret this to mean that your procedure might be (I ain't betting on it) a good approximation for $\rho_{00}$, $\rho_{01}$, and $\rho_{11}$, but it's a terrible approximation on any matrix element with a 2 in it.

If you want to check this, try incrementally raising the dimension on the matrix and seeing if $\rho_{00}$ converges. Hopefully, the negative eigenvalue will converge to a non-negative eigenvalue. If not, it might be worth investigating L more carefully.

Of course. Then you deal with an "open quantum system". One approach is the Lindblad equation. See, e.g., Weinberg, Lectures on Quantum Mechanics.

 

[Haorong Wu]

Hi, @Twigg . I would like to share my latest result. I have enlarged the dimension to five, including $\left | 3 \right >$ and $\left | 4 \right >$.

The result shows that there are two negative eigenvalues, -0.06 and -0.01, now. The states associated with negative eigenvalues have high coefficients before $\left | 4 \right >$.

Hence, I think your statement is correct. My model is not a well approximation for the highest states. For a state $\left | p \right >$, it directly couples to $\left | p-1 \right >$, and $\left | p+1 \right >$. However, for a highest state in the subspace, only $\left |p-1 \right >$ is considered, so half of coupling is missing.

But given that the negative eigenvalues are quite small compared to the largest positive eigenvalue, 0.9622, I believe that it is safe to throw those states associated with negative eigenvalues away.

 

[A. Neumaier]

@A. Neumaier When I first studied decoherence, I tried to evolve a state and then calculate the density matrix. However, when a state $\left | \psi \right >$ evolves into another state$\left | \psi ' \right >$, the density matrix $\left | \psi ' \right >\left < \psi ' \right |$ will always be a pure state, as @Twigg pointed out. No decoherence can be deduced from them.

@A. Neumaier I believe the OP is saying that the system is expected to evolve into a mixed state and $\psi$ is not defined. Am I misunderstanding something here?

In my reply at #7, $\psi$ is a matrix representing the mixed state, not a pure vector state!

 

[Haorong Wu]

Hi, @A. Neumaier . I am sorry I am not familiar with the representation of mixed states by matrix. Do you konw any textbook about this? Thanks!

 

[A. Neumaier]

This is simple linear algebra: Any positive semidefinite Hermitan matrix $A$ can be written as $A=LL^*$; e.g., by a Cholesky factorization ($L$ triangular) or a spectral factorization ($L$ Hermitian positive semidefinite). Doing this with $\rho$ in place of $A$ and $\psi$ in place of $L$, you can check that the differential equation I wrote for $\psi$ leads to you equation for $\rho$. So just factor your initial $\rho$ to get the initial $\psi$, then solve the 9 equations for the matrix $\psi$ using ode45 with 9 variables, and calculate from the $\psi(t)$ the $\rho(t)$!

 

[Haorong Wu]

Thanks! @A. Neumaier . I will try that method now. I am sorry I misunderstood your idea. When I saw $\psi$, the state vectors immediately pop out in my head. lol

 

[A. Neumaier]

@A. Neumaier when a state $\left | \psi \right >$ evolves into another state$\left | \psi ' \right >$, the density matrix $\left | \psi ' \right >\left < \psi ' \right |$ will always be a pure state, as @Twigg pointed out.

With your differential equation, density matrices representing pure states will in fact remain pure. So there is no decoherence because of your assumed dynamics!

To get decoherence you need to add double commutator terms to your dynamics for $\rho$. See, e.g. https://en.wikipedia.org/wiki/Lindbladian

 

[f95toli]

I am definitely not an expert, but to me it sounds like the OPs approvach is more similar to Bloch-Redfield than Lindblad, or?

 

[A. Neumaier]

With your differential equation, density matrices representing pure states will in fact remain pure. So there is no decoherence because of your assumed dynamics!

In fact, your initial density operator belongs to a pure state, hence with an exact integration you should get a rank 1 density matrix at all times, i.e., the two small eigenvalues are spurious, due to discretization errors. Decrease the step size of the ODE solvers to verify this!

 

[Haorong Wu]

With your differential equation, density matrices representing pure states will in fact remain pure. So there is no decoherence because of your assumed dynamics!

To get decoherence you need to add double commutator terms to your dynamics for $\rho$. See, e.g. https://en.wikipedia.org/wiki/Lindbladian

I am not sure how to put it correctly. Well, what I am studying is not pure quantum-mechanics. Maybe semi-QM. And $L_{mn}$ has nothing to do with the Hamiltonian, so the evolution is not unitary. Instead, it has a quite ugly form, e.g.,
$$L_{mn}=\frac {m\sqrt{(|l|+r+1)(r+1)}(kw_0^2+2iz)}{8k^2w_z^2\tilde z^4}[ 6iz+8k\tilde z(z+\tilde z)+kw_0^2(3-4ik\tilde z)] $$ for $p_m=r+1$ and $p_n=r$.
So decoherence may happen in this case.

 

[Haorong Wu]

I am definitely not an expert, but to me it sounds like the OPs approvach is more similar to Bloch-Redfield than Lindblad, or?

Thanks! I think I will get to know them.

 

[A. Neumaier]

I am not sure how to put it correctly. Well, what I am studying is not pure quantum-mechanics. Maybe semi-QM. And $L_{mn}$ has nothing to do with the Hamiltonian, so the evolution is not unitary. Instead, it has a quite ugly form, e.g.,
$$L_{mn}=\frac {m\sqrt{(|l|+r+1)(r+1)}(kw_0^2+2iz)}{8k^2w_z^2\tilde z^4}[ 6iz+8k\tilde z(z+\tilde z)+kw_0^2(3-4ik\tilde z)] $$ for $p_m=r+1$ and $p_n=r$.
So decoherence may happen in this case.

What you have is an effective, nonhermitian Hamiltonian, thus there is dissipation (the trace of the density matrix decreases). But there is no decoherence, independent of the form of $L$. To see this, start in my equation with a vector $\psi$ rather than a matrix, and verify that the associated density matrix of rank 1 satisfies your equation. Since the initial value problem has a unique solution, the rank 1 solution is the only one. Thus pure initial states remain pure at all times.

Indeed, effective Hamiltonians are characteristic for simple dissipative problems in which decoherence can be neglected.

 

[Haorong Wu]

What you have is an effective, nonhermitian Hamiltonian, characteristic for simple dissipative problems. Thus there is dissipation (the trace of the density matrix decreases) But there is no decoherence, independent of the form of $L$. To see this, start in my equation with a vector $\psi$ rather than a matrix, and verify that the associated density matrix of rank 1 satisfies your equation. Since the initial value problem has a unique solution, the rank 1 solution is the only one. Thus pure initial states remain pure at all times.

Well, I have tried it in a slightly different scenario. Strange thing happened. The trace of $\rho=\psi \psi^*$ decreases to almost zero quite fast. However, if I use the matrix to calculate, the trace of $\rho$ remains almost 1.

Something must be wrong. I have some suspects, and I will dig it further tonight.

 

[Twigg]

@A. Neumaier Sorry, I was on a road trip and wasn't able to reply. Thanks for explaining what you meant about $\psi$ being a matrix and the equation being dissipative but not decoherent. There were more subtleties here than I realized! Much appreciated