如何在 python 上解决 TISE 的简单边值问题

我正在尝试求解无限势阱的 TISEV=0在间隔上[0,L]。这个练习给我们提供了波函数及其导数的值0 is 0,1分别。这使我们能够使用scipy.integrate.odeint函数来解决给定能量值的问题。

现在的任务是在给定波函数的进一步边界条件下找到能量特征值L is 0，使用 python 上的求根函数。我做了一些研究，只能找到一种叫做“拍摄方法”的东西，但我不知道如何实施。另外，我遇​​到了求解 BVP scipy 函数，但是我似乎无法理解该函数的第二个输入到底是什么（边界条件残差）

m_el   = 9.1094e-31      # mass of electron in [kg]
hbar   = 1.0546e-34      # Planck's constant over 2 pi [Js]
e_el   = 1.6022e-19      # electron charge in [C]
L_bohr = 5.2918e-11      # Bohr radius [m]

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

def eqn(y, x, energy):              #array of first order ODE's     
    y0 = y[1]
    y1 = -2*m_el*energy*y[0]/hbar**2

    return np.array([y0,y1])

def solve(energy, func):           #use of odeint
    p0 = 0
    dp0 = 1
    x = np.linspace(0,L_bohr,1000)
    init = np.array([p0,dp0])
    ysolve = odeint(func, init, x, args=(energy,))
    return ysolve[-1,0]

这里的方法是在solve(energy,func)中输入eqn作为func。 L_bohr是这个问题中的L值。我们正在尝试使用某种 scipy 方法以数值方式找到能量特征值

对于 scipy 中的所有其他求解器，参数顺序x,y，甚至在odeint可以通过给出选项来使用此命令tfirst=True。从而更改为

def eqn(x, y, energy):              #array of first order ODE's     
    y0, y1 = y
    y2 = -2*m_el*energy*y0/hbar**2

    return [y1,y2]

对于 BVP 求解器，您必须将能量参数视为 具有零导数的额外状态分量，从而添加第三个槽 在边界条件中。西皮的solve_bvp允许将其保留为参数， 这样你就可以在边界条件中得到 3 个槽，从而可以将一阶导数固定在x=0从特征空间中选择一个非平凡的解决方案。

def bc(y0, yL, E):
    return [ y0[0], y0[1]-1, yL[0] ]

接下来构造一个接近疑似基态的初始状态并调用求解器

x0 = np.linspace(0,L_bohr,6);
y0 = [ x0*(1-x0/L_bohr), 1-2*x0/L_bohr ]
E0 = 134*e_el

sol = solve_bvp(eqn, bc, x0, y0, p=[E0])
print(sol.message, "  E=", sol.p[0]/e_el," eV")

然后产生情节

x = np.linspace(0,L_bohr,1000)
plt.plot(x/L_bohr, sol.sol(x)[0]/L_bohr,'-+', ms=1)
plt.grid()

The algorithm converged to the desired accuracy. E= 134.29310361903723 eV