我是python的新手，正在为我的物理项目做编码，这需要生成具有变量E的矩阵

请注意以后的问题，请不要使用诸如numpy之类的未使用进口内容，而不要像# a = eye(3,3)之类的僵尸代码。这有助于保持代码尽可能简洁和简短。另外，由于缩进问题，示例代码将无法运行，因此在复制和粘贴代码时，请先确保它可以工作。始终尝试使您的问题尽可能简短和模块化。

T[0,0]的表达式过于复杂，无法用SymPy进行解析求解，因此需要数值逼近。剩下2个选项：

1. 使用先进的SciPy求解器，但由于SciPy不会以任何方式处理SymPy对象，因此需要将类型强制转换为float值。
2. 使用SymPy的根求解器，该求解器虽然不太先进，但使用起来可能更简单。

这两种方法都只会产生一个数字作为输出，因为您不能期望数值求解器会找到每个根。如果要查找多个点，则建议您使用要用作初始值的点列表，将每个点输入到求解器中并跟踪不同的输出。但是，这绝不能保证您已获得每个根。

只有在没有问题的情况下使用SciPy和SymPy才可以混合使用。 SciPy根本无法与SymPy一起玩，在使用SciPy时，您应该只具有list，float和complex实例。

import math
import sympy as sy
from scipy.optimize import newton

# Constants(Values at temp 10K)
hbar = 1.055E-34
m0 = 9.1095E-31  # free mass of electron
q = 1.602E-19
v = [0.510, 0, 0.510]  # conduction band offset in eV
m1 = 0.043  # effective mass in In_0.53Ga_0.47As
m2 = 0.072  # effective mass in Al_0.48In_0.52As
d = [-math.inf, 100, math.inf]  # dimension of structure in nanometers

'''scaling factor to with units of E in eV, mass in terms of free mass of electron, length in terms 
of nanometers '''
s = (2 * q * m0 * 1E-18) / hbar ** 2
E = sy.symbols('E')  # Suppose energy of incoming particle is 0.3eV
m = [0.043, 0.072, 0.043]  # effective mass of electrons in layers
for i in range(3):
    print('Effective mass of e in layer', i, 'is', m[i])
k = []  # Defining an array for wavevectors in different layers
for i in range(3):
    k.append(sy.sqrt(s * m[i] * (E - v[i])))
    print('Wave vector in layer', i, 'is', k[i])

x = []
for i in range(2):
    x.append((k[i + 1] * m[i]) / (k[i] * m[i + 1]))
    # Define Boundary condition matrix for two interfaces.
D0 = (1 / 2) * sy.Matrix([[1 + x[0], 1 - x[0]], [1 - x[0], 1 + x[0]]], dtype=complex)
D1 = (1 / 2) * sy.Matrix([[1 + x[1], 1 - x[1]], [1 - x[1], 1 + x[1]]], dtype=complex)

# Define Propagation matrix for 2nd layer or quantum well
P1 = 1 * sy.Matrix([[sy.exp(-1j * k[1] * d[1]), 0], [0, sy.exp(1j * k[1] * d[1])]], dtype=complex)
print("abs")
T = D0 * P1 * D1

# did not converge for 0.5
x0 = 0.75

# method 1:
def f(e):
    # evaluate T[0,0] at e and remove all sympy related things.
    result = complex(T[0, 0].replace(E, e))
    return result

solution1 = newton(f, x0)
print(solution1)

# method 2:
solution2 = sy.nsolve(T[0,0], E, x0)
print(solution2)

此打印：

(0.7533104353644469-0.023775286117722193j)
1.00808496181754 - 0.0444042144405285*I

请注意，第一行是本地Python complex实例，而第二行是SymPy复数的实例。一个人可以简单地用print(complex(solution2))来转换第二个。

现在，您会注意到它们产生的数字不同，但是两者都是正确的。从Geogebra图可以看出，此函数似乎有很多零：

红色轴为Re(E)，绿色轴为Im(E)，蓝色轴为|T[0,0]|。这些“峰值”中的每一个可能都是零。