在 sympy 中替换一个矩阵

Question

我正在使用 pythTB 研究石墨烯的紧密结合模型。我想在计算中加入自旋元素。 rashba 跳跃项的 hamiltonian 使 pauli 自旋矩阵向量与站点跳跃向量交叉。

最初我创建了一个矩阵列表并将其与向量相交，不幸的是，这并没有产生正确的结果（我认为在取向量叉积之后，再取矩阵的叉积）。

接下来，我声明了 3 个符号“s_x”、“s_y”和“s_z”，并使用它们代替了我的 pauli 自旋矩阵向量中的矩阵。服用叉积后，我收到了正确的结果。我遇到的问题是我无法将矩阵替换为我添加的变量符号。可以这样做吗？还是我需要手动取叉积？

这是我的一些代码：

from __future__ import print_function
from pythtb import * # import TB model class
from sympy import symbols
import numpy as np
import matplotlib.pyplot as plt

# create list of pauli spin matrices 
sx = [[0., 1.],[1., 0.]]
sy = [[0., -1.j],[1.j, 0.]]
sz = [[1., 0.],[0., -1.]]
Id = [[1., 0.], [0., 1.]]
s_pauli = np.zeros((4, 2, 2), dtype=complex)
s_pauli = [Id, sx, sy, sz]

# create s_pauli without identity matrix
s_pau = np.zeros((3, 2, 2), dtype=complex)
s_pau = [ s_x, s_y, s_z]

ab00 = [ 0.5, 0.28867513, 0.]

sig_x_ab00 = np.cross( s_pau, ab00)

如果我打印sig_x_ab00[2]（这是我目前唯一感兴趣的），那么我会得到：

0.288675134594813*s_x - 0.5*s_y

获得该信息后，我想通过执行以下命令将s_pauli[1] 替换为s_x 和s_pauli[2] 替换s_y：

sig_x_ab00_ = sig_x_ab00.subs(s_x, s_pauli[1])

我得到以下错误输出：

AttributeError: 'numpy.ndarray' object has no attribute 'subs'

我所做的一切都有效吗？还是有更好的方法来解决这个问题？

非常感谢任何输入！ 谢谢！

Answer 1

让我们运行您的代码，但要查看每一步。不要做任何假设。

我正在使用isympy 交互环境； ipython 与 sympy 增强。我还导入了np。

In [4]: ab00 = [ 0.5, 0.28867513, 0.]                                           

In [5]: s_pauli                                                                 
Out[5]: 
[[[1.0, 0.0], [0.0, 1.0]],
 [[0.0, 1.0], [1.0, 0.0]],
 [[0.0, (-0-1j)], [1j, 0.0]],
 [[1.0, 0.0], [0.0, -1.0]]]

这是一个列表。前面的np.zeros(...) 表达式什么也不做。在 Python 中，我们不设置变量的“类型”。

我们可以从这个列表中创建一个数组：

In [6]: np.array(s_pauli)                                                       

s_pauli[1] 有效，因为它只是列表索引。

以及添加的符号：

In [11]: s_x, s_y, s_z = symbols('s_x s_y s_z')                                 

In [12]: s_x                                                                    
Out[12]: sₓ

In [13]: s_pau = [ s_x, s_y, s_z]                                               

同样，s_pau 是一个列表，而不是一个数组。在cross中使用时会变成一个数组：

In [14]: np.array(s_pau)                                                        
Out[14]: array([s_x, s_y, s_z], dtype=object)

请注意，这是一个对象 dtype 数组，它仍然非常像一个列表。一些基本的数学工作，因为像乘法和加法这样的数学是为符号定义的。但是像 np.log 和 np.sin 这样的先验函数在这样的数组上不起作用。

cross 只使用乘法和加法，所以它适用于这些对象数组：

In [15]: sig = np.cross( s_pau, ab00)                                           

In [16]: sig                                                                    
Out[16]: array([-0.28867513*s_z, 0.5*s_z, 0.28867513*s_x - 0.5*s_y], dtype=object)

sig 是一个 numpy 数组。它不是一个 sympy 表达式，也没有 subs 方法。同样，密切关注正在发生的事情是值得的。

数组的元素是 sympy 表达式：

In [17]: sig[2]                                                                 
Out[17]: 0.28867513⋅sₓ - 0.5⋅s_y

In [20]: s2 = sig[2]                                                            

具有标量值的 subs 有效：

In [22]: s2.subs(s_x, 1)                                                        
Out[22]: 0.28867513 - 0.5⋅s_y

但不是列表

In [23]: s2.subs(s_x, s_pauli[1])                                               
Out[23]: 0.28867513⋅sₓ - 0.5⋅s_y

但是，如果我从中制作 sympy 矩阵：

In [24]: s_pauli[1]                                                             
Out[24]: [[0.0, 1.0], [1.0, 0.0]]

In [25]: Matrix(s_pauli[1])                                                     
Out[25]: 
⎡0.0  1.0⎤
⎢        ⎥
⎣1.0  0.0⎦

In [26]: s2.subs(s_x, Out[25])                                                  
Out[26]: 
           ⎡    0       0.28867513⎤
-0.5⋅s_y + ⎢                      ⎥
           ⎣0.28867513      0     ⎦

替换确实有效。

通常混合sympy 和numpy 是命中注定的；一些工作，几乎更多的是偶然而不是设计。其他人没有。 sympy.lambdify 是制作可与 numpy 数组一起使用的函数的最可靠方法。

在这种情况下，我怀疑您最好使用 cross 的 sympy 版本，并进行 sympy.Matrix 替换。