Logarithme des matrices singulières
Je voulais calculer la quantité suivante:
$X=Tr[\rho_1 \log[\rho_2]]$,
comme dans l' entropie relative . Ici,$\rho_1, \rho_2$ sont des matrices semi-définies positives avec support non orthogonal (pour que la chose ne diverge pas - des termes comme $0 \log[0]$ devrait être ramené à zéro, ce qui est une hypothèse standard) et $Tr$ est trace.
Il y a une question SE similaire , la fonction donnée
MatrixLogSafe[x_] := MatrixFunction[Piecewise[{{Log[#1], #1 > 0}}] &, x]
qui devrait traiter du logarithme de la matrice, il se comporte cependant de façon étrange.
Par ex je suppose que $\rho_1 =\rho_2 =\{\{0.33,0,0\},\{0,0,0\},\{0,0,0.66\}\}$. La quantité$X$ devrait alors être
$X=0.33 \log[0.33] + 0.66 \log[0.66]= -0.640099$.
Cependant, l'utilisation de MatriLogSafe dans la définition donne une sortie différente:
In[402]:= Tr[{{0.33, 0, 0}, {0, 0, 0}, {0, 0, 0.66}}.MatrixLogSafe[{{0.33, 0, 0}, {0, 0, 0}, {0, 0, 0.66}}]]
Out[402]= -0.731717
Le problème est que MatrixLogSafe "change parfois les vecteurs propres",
In[403]:= MatrixLogSafe[{{0.33, 0, 0}, {0, 0, 0}, {0, 0, 0.66}}]
Out[403]= {{0., 0., 0.}, {0., -0.415515, 0.}, {0., 0., -1.10866}}
(alors $\log[0.33]= -1.10866$ et $\log[0.66]=-0.415515$, mais la sortie doit être {{-1.10866, 0, 0.}, {0., 0., 0.}, {0., 0., - 0.415515}}).
(D'une certaine façon, je pense que le problème est que j'utilise des valeurs numériques, mais je veux que la fonction fonctionne à la fois pour les nombres numériques et "exacts" (?))
Comment y remédier?
J'ai considéré la réponse donnée par Carl Woll, mais quelque chose ne fonctionne toujours pas. En particulier, considérons deux matrices,$\rho_1$:
{{1/4, 1/4 E^(-((I \[Pi])/10)), 1/4 E^(-((I \[Pi])/10)), 1/
4}, {1/4 E^((I \[Pi])/10), 1/4, 1/4,
1/4 E^((I \[Pi])/10)}, {1/4 E^((I \[Pi])/10), 1/4, 1/4,
1/4 E^((I \[Pi])/10)}, {1/4, 1/4 E^(-((I \[Pi])/10)),
1/4 E^(-((I \[Pi])/10)), 1/4}}
et $\rho_2$:
{{1/4 (1 - 1/Sqrt[E]) + 1/(4 Sqrt[E]),
1/4 E^(-(1/2) - (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10)),
1/4 E^(-(1/2) - (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10)),
1/(4 Sqrt[
E]) + (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10))^2}, {1/4 E^(-(1/2) + (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10)), 1/4 (1 - 1/Sqrt[E]) + 1/(4 Sqrt[E]),
1/(4 Sqrt[
E]) + (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10))^2,
1/4 E^(-(1/2) + (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10))}, {1/4 E^(-(1/2) + (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10)),
1/(4 Sqrt[
E]) + (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10))^2, 1/4 (1 - 1/Sqrt[E]) + 1/(4 Sqrt[E]),
1/4 E^(-(1/2) + (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10))}, {1/(
4 Sqrt[E]) + (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10))^2,
1/4 E^(-(1/2) - (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10)),
1/4 E^(-(1/2) - (I \[Pi])/10) +
1/2 (1 - 1/Sqrt[E]) (1/4 E^(-((I \[Pi])/10)) +
1/4 E^((I \[Pi])/10)), 1/4 (1 - 1/Sqrt[E]) + 1/(4 Sqrt[E])}}
qui sont peut-être laids, mais ils sont tous les deux hermitiens et positifs-semi-définis.
Je deviens cependant complexe $X$ ($X=densityTrace[\rho_1,\rho_2]=-0.0019613 + 0.393667 I$), ce qui ne peut pas être!
Cependant, si je calcule d'abord les matrices numériques, j'obtiens une sortie réelle $densityTrace[N[\rho_1],N[\rho_2]]=-0.0432473$ (?!)
J'ai légèrement changé la définition de la densitéTrace
densityTrace[a_, b_] := Module[{λ, S, d},
{λ, S} = Eigensystem[b];
S = Transpose[S];
d = Diagonal[Inverse[S]. a. S];
Total @ MapThread[If[Chop[#1]==0,0,Chop[#1] Log[Chop[#2]]]&, {d, λ}]
]
(ajoutant Chop), pour se débarrasser de quelques très petits «déchets» imaginaires, est-ce une bonne idée?
Réponses
Peut-être que vous pouvez utiliser Eigensystempour transformer votre matrice singulière$\rho_2$ dans $S.d.S^{-1}$, puis utilisez la propriété cyclique de Tr. C'est:
$$ \begin{align} \text{Tr}{\left(\rho_1 . \log(\rho_2) \right)}&= \text{Tr}{\left(\rho_1 . \log(S . d . S^{-1})\right)} \\ &= \text{Tr}{\left(\rho_1 . S . \log(d) . S^{-1} \right)}\\ &= \text{Tr}{\left(S^{-1} . \rho_1 . S . \log(d)\right)}\\ &= \text{Diagonal}{\left( S^{-1} . \rho_1 . S\right) . \log (d)} \end{align} $$
Un code qui implémente cette idée:
densityTrace[a_, b_] := Module[{λ, S, d},
{λ, S} = Eigensystem[b];
S = Transpose[S];
d = Diagonal[Inverse[S]. a. S];
Total @ MapThread[If[#1==0,0,#1 Log[#2]]&, {d, λ}]
]
Pour vos exemples de matrices:
m1 = m2 = DiagonalMatrix[{.33, 0, .66}];
densityTrace[m1, m2]
-0,640099
Il est probablement plus facile de détecter la singularité plus tôt:
MatrixLogSafe[x_?SquareMatrixQ] /; Det[x] > 0 := MatrixLog[x];
MatrixLogSafe[x_?SquareMatrixQ] := ConstantArray[0, Dimensions[x]];
Je peux me tromper sur la condition exacte dont vous avez besoin. En tout cas, je ne m'attends pas MatrixFunction[Piecewise[...]]à vous faire beaucoup de bien. Je ne sais même pas ce que la condition # > 0est censée faire à l'intérieur MatrixFunction.