2
$\begingroup$

I am trying to solve the following Lindblad equation analytically with Mathematica:

$$\frac{d\rho}{dt}=-i[H,\rho]+\sum\nolimits_{i=0}^1 K_i\rho K_i^+-\frac12\cdot\{K_i^+K_i,\rho\},$$ where $$H=\begin{pmatrix} -1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix},$$

$$K_1=\frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix},$$

$$K_2=\frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 0 & 0 & 0\\ 0 & -1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & -1 \end{pmatrix}.$$

My code looks like this:

Clear["Global`*"]
H[t_] = {{-1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, 1, 0}, {0, 0, 0, -1}};
K1 = {{1/Sqrt[2], 0, 0, 0}, {0, 1/Sqrt[2], 0, 0}, {0, 0, -1/Sqrt[2], 
    0}, {0, 0, 0, -1/Sqrt[2]}};
K1plus = ConjugateTranspose[K1];
K2 = {{1/Sqrt[2], 0, 0, 0}, {0, -1/Sqrt[2], 0, 0}, {0, 0, 1/Sqrt[2], 
    0}, {0, 0, 0, -1/Sqrt[2]}};
K2plus = ConjugateTranspose[K2];
eqn = I*D[rho[t], t] == -I*H[t] . rho[t] + I*rho[t] . H[t] + 
    K1 . rho[t] . K1plus - 1/2*K1plus . K1 . rho[t] - 
    1/2*rho[t] . K1plus . K1 + K2 . rho[t] . K2plus - 
    1/2*K2plus . K2 . rho[t] - 1/2*rho[t] . K2plus . K2;

sol = DSolve[eqn, rho, t]

Nevertheless, I don't get the exact solution. Mathematica only prints the input DSolve[...].

Does that mean that it can't find the solution?

$\endgroup$
0

1 Answer 1

3
$\begingroup$

DSolve can't automatically determine the dimensionality of your solution. Instead it assumes it to be scalar. To solve it as a matrix, you have to provide it as a matrix $$\rho(t)=\begin{pmatrix}\rho_{1,1}(t) & \rho_{1,2}(t) & \dots\\\rho_{2,1}(t) & \rho_{2,2}(t) & \dots \\\vdots & \vdots & \ddots\end{pmatrix}$$

You can do this easily with Table:

rhoMat = Table[rho[i][j][t], {i, 1, Length@K1}, {j, 1, Length@K1}];

The replace every occurence of rho[t] with rhoMat:

eqn = I*D[rhoMat, t] == -I*H[t] . rhoMat + I*rhoMat . H[t] + 
    K1 . rhoMat . K1plus - 1/2*K1plus . K1 . rhoMat - 
    1/2*rhoMat . K1plus . K1 + K2 . rhoMat . K2plus - 
    1/2*K2plus . K2 . rhoMat - 1/2*rhoMat . K2plus . K2;

Plug this into DSolve and get your solutions:

DSolve[eqn // Flatten, rhoMat // Flatten, t];
First[rhoMat /. DSolve[eqn // Flatten, rhoMat // Flatten, t]]

{{C[1], E^((2 + I) t) C[2], E^((2 + I) t) C[3], E^(2 I t) C[4]}, {E^((-2 + I) t) C[5], C[6], E^(2 I t) C[7], E^((-2 + I) t) C[8]}, {E^((-2 + I) t) C[9], E^(2 I t) C[10], C[11], E^((-2 + I) t) C[12]}, {E^(2 I t) C[13], E^((2 + I) t) C[14], E^((2 + I) t) C[15], C[16]}}

Btw, I think you added an unnecessary I infront of the time derivative. The original Lindbladian does not feature this.


EDIT: When I solve the Lindbladian, I usually do not take DSolve. Instead I reduce the problem to $$\rho'(t)=M\rho(t)$$ Which can easily be done with CoefficientArrays and then solved as a Eigensystem problem. I paste the code for your equation here (ommiting the I whih I think is wrong).

eqn = -I*H[t] . rhoMat + I*rhoMat . H[t] + K1 . rhoMat . K1plus - 
   1/2*K1plus . K1 . rhoMat - 1/2*rhoMat . K1plus . K1 + 
   K2 . rhoMat . K2plus - 1/2*K2plus . K2 . rhoMat - 
   1/2*rhoMat . K2plus . K2;

systemMatrix = Last[CoefficientArrays[Flatten@eqn, Flatten@rhoMat]];

{eigVec, eigSys} = Eigensystem[mat2];
sol = Total[
   Table[(Exp[eigVec[[i]]*t]*eigSys[[i]]*C[i]), {i, 1, 
     Length@eigVec}]];
sol = rhoMat /. Thread[Flatten[rhoMat] -> sol]

{{C[16], E^((-1 + 2 I) t) C[4], E^((-1 + 2 I) t) C[3], E^(-2 t) C[12]}, {E^((-1 - 2 I) t) C[8], C[15], E^(-2 t) C[11], E^((-1 - 2 I) t) C[7]}, {E^((-1 - 2 I) t) C[6], E^(-2 t) C[10], C[ 14], E^((-1

  • 2 I) t) C[5]}, {E^(-2 t) C[9], E^((-1 + 2 I) t) C[2], E^((-1 + 2 I) t) C[1], C[13]}}

Besides the enumeration of the constants, the solutions are identical. With this approach, one can also easily solve for initial value problems with $$\rho(0)=m_{eig}^TC$$

For which $C$ is your sought after vector of constants, $m_{eig}$ the eigensystem matrix and $\rho(0)$ your intial boundary values.

$\endgroup$

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.