Matplotlib - 绘制具有相同尺寸的 2 个值时出错

Question

def Hamiltonian(alpha,h):

    Sx = np.array([[0,1],[1,0]])
    Sy = np.array([[0,-1j],[1j,0]])
    Sz = np.array([[1,0],[0,-1]])
    I  = np.array([[1,0],[0,1]])

    H = -1*((alpha*np.kron(np.kron(Sx,Sx),I))
       + (alpha*np.kron(np.kron(Sy,Sy),I))
       + (alpha*np.kron(np.kron(I,Sx),Sx))
       + (alpha*np.kron(np.kron(I,Sy),Sy))
       + (h*np.kron(np.kron(I,Sz),I)))

    return H
np.set_printoptions(linewidth=100)
Hamiltonian(1,0.5).real

它返回以下矩阵（为了清楚起见，我在下面尝试做的只是输入这个）

定义哈密顿量后，我想寻找熵作为 h 参数的函数。对于这类问题，代码背后的物理原理并不重要。

# von Neumann entropy as a function of h and beta - Complete
# Definition of a mixed state: [Thermal Density Matrix used]
h = np.arange(0,2.5,0.1)
beta = 2
for i in range(h.size):
    H = Hamiltonian(1.0, h[i] )
    rho_thermal = expm(-1.0 * beta * H)
    tr = np.trace(rho_thermal)
    rho_thermal = rho_thermal / tr
    np.set_printoptions(linewidth=100)
    eigvals_rho_thermal, eigvecs_rho_thermal = LA.eigh(rho_thermal)
    # Entropy
    s = 0.0
    for i in range(eigvals_rho_thermal.size):
        s += -1.0 * (eigvals_rho_thermal[i] * np.log(eigvals_rho_thermal[i]))
    print(s)
plt.plot(h,s)

我的问题是为什么会出现以下错误：

我的代码为 s 返回 25 个值，为 h 返回 25 个值，为什么不绘制它们？

Answer 1

您需要将在每个循环步骤中获得的s 的值存储在某处，以便以后能够为每个h 绘制每个s。通常，您可以稍微简化代码并使用函数而不是 for 循环。 以下是我的做法。

import numpy as np
import matplotlib.pyplot as plt
from scipy.linalg import expm

def Hamiltonian(alpha,h):

    Sx = np.array([[0,1],[1,0]])
    Sy = np.array([[0,-1j],[1j,0]])
    Sz = np.array([[1,0],[0,-1]])
    I  = np.array([[1,0],[0,1]])

    H = -1*((alpha*np.kron(np.kron(Sx,Sx),I))
       + (alpha*np.kron(np.kron(Sy,Sy),I))
       + (alpha*np.kron(np.kron(I,Sx),Sx))
       + (alpha*np.kron(np.kron(I,Sy),Sy))
       + (h*np.kron(np.kron(I,Sz),I)))

    return H

# von Neumann entropy as a function of h and beta - Complete
# Definition of a mixed state: [Thermal Density Matrix used]
def get_entropy(beta, h, alpha=1.0):
    H = Hamiltonian(alpha, h)
    rho_thermal = expm(-1.0 * beta * H)
    tr = np.trace(rho_thermal)
    rho_thermal = rho_thermal / tr
    eigvals_rho_thermal, eigvecs_rho_thermal = np.linalg.eigh(rho_thermal)
    # Entropy
    s = -np.sum(eigvals_rho_thermal*np.log(eigvals_rho_thermal))
    return s

h = np.arange(0,2.5,0.1)
beta = 2    
s = [get_entropy(beta, hi) for hi in h]
plt.plot(h,s)

plt.show()