量子谐振子/简谐振子模型 [英] Modelling a quantum harmonic oscillator/SHM

问题描述

我需要帮助来找出为什么我的b)的基本状态图看起来是错误的，这里是完整的问题：(我认为完整地发布它会给我试图使用的方法提供上下文)

(A)考虑两个无限高的壁之间的正方形势垒，𝑉(𝑥)=0，距离𝐿等于玻尔半径，即对于区间[0，𝐿]中的所有x。

• 编写一个接受参数Energy的函数solve(energy, func)和一个python函数𝑓𝑢𝑛𝑐。此函数应解决上述情况下的薛定谔模式，并仅返回边界𝜓(𝐿)在边界𝐿的最终值。

• 使用函数solve(energy, func)编写脚本，以EV为单位计算该势垒中电子的基态能量(即，将结果除以基本电荷值)。有关初始条件，请参阅下面的技术提示。对于要为𝜓(𝑥求解的值数)，建议使用值1,000。

您的计算结果应该是一个介于134 eV和135 eV之间的数字。其中一项测试将评估您的解决方案(能量，函数)函数的扭曲势井。

(B)考虑调和势 𝑉(𝑥)=𝑉₀𝑥²/𝑎²，，其中𝑉₀和𝑎=10⁻¹¹m是常量。将变量𝑥的无限范围限制为区间[−10𝑎，10𝑎]，𝑉₀=50 eV。

• 谐振子有等距的能量本征值。通过计算基态和前两个激发态来检查这一点是否属实，以达到计算的精度。(提示：基态的能量范围在100到200 eV之间。)

  /li>

• 为了测试您的结果，请编写一个函数result()，该函数必须将EV中计算出的能量本征值的差值作为单个数字返回。请注意，具有预期数量的测试是隐藏的，并且将以±1 eV的精度测试您的结果。

• 提供绘图标题和相应的轴标签，使用彩色编码线将所有三个波函数绘制到单个画布上，并提供适当的x和y轴限制以使所有曲线清晰可见。

技术提示：这不是薛定谔颂歌的初值问题，而是边值问题！这需要付出额外的努力，采用与前面的颂歌练习不同的方法，因此是一个"看不见"的问题。

• 分别为𝑥₀=0或𝑥₀=−10𝑎的两个问题设置一个简单的初始条件：𝜓(𝑥₀)=0和𝑑𝜓(𝑥₀)/𝑑𝑥=1。以此为起点求解常微分方程组，并改变能量𝐸，直到解收敛。任务是评估能量变量的变化，直到满足第二边界条件(例如，练习(A)的L)，因为由于初始条件的选择，第一边界条件已经满足。

• 这需要对解可能所在的能量区间的初始猜测和求根的计算方法。在Scipy中搜索求根方法，并选择一种不需要导数知识的方法。这里适合根查找，因为要满足的边界条件是𝜓(𝑥)=0。

量子物理背景，这两个练习的边界条件都是在每个势边界上𝜓(𝑥=0。这些仅在特定的离散能量值𝐸下实现，称为能量本征值，其中最小的称为基态能量。

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]
V0 = 50*e_el
a = 10**(-11)

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

def eqn(y, x, energy):                       #part a)
    y0 = y[1]
    y1 = -2*m_el*energy*y[0]/hbar**2

    return np.array([y0,y1])

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

def eigen(energy):

    return solve(energy, eqn)

root_ = optimize.toms748(eigen,134*e_el,135*e_el)
root = root_/e_el

print('Ground state infinite square well',root,'eV')

intervalb = np.linspace(-10*a,10*a,1000)                    #part b)

def heqn(y, x2, energy):
    f0 = y[1]
    f1 = (2.0 * m_el / hbar**2) * (V0 * (x2**2/a**2) - energy) * y[0]

    return np.array([f0,f1])

def solveh(energy, func):
    ph0 = 0
    dph = 1
    init = np.array([ph0,dph])
    ysolve = odeint(func, init, intervalb, args=(energy,))
    return ysolve

def boundary(energy):        #finding the boundary V=E to apply the b.c
    f = a*np.sqrt(energy/V0)
    index = np.argmin(np.abs(intervalb-f))

    return index

def eigen2(energy):

    return solveh(energy,heqn)[boundary(energy),0]

groundh_ = optimize.toms748(eigen2,100*e_el,200*e_el)
groundh = groundh_/e_el


print('Ground state of Harmonic Potential:', groundh, 'eV')
plt.suptitle('Harmonic Potential Well')
plt.xlabel('x (a [pm])')
plt.ylabel('Psi(x)')

groundsol = solveh(groundh_,heqn)[:,0]
plt.plot(intervalb/a, groundsol)

对于100 eV到200 eV之间的所有能量值，图表形状看起来都是这样的。我不知道我做错了什么。我已经尽可能地测试了我的代码。

推荐答案

函数的代码或任务文本没有任何原因boundary。如a)取边界条件的正确区间末端的值，理想情况下这应该是无穷大的值。

代码中还有另外两个问题，它们实际上不是您的错误：

• 由于scipy的最新版本中可能没有的某些原因，并且我在可用的源代码中找不到它，您使用的根查找器toms748的调用似乎缓存了函数，而不是用新的函数替换它。这意味着在部分b)中，根查找器调用仍然找到根，但它是来自部分a)的根。只需检查函数值即可确认这一点。我建议对初始括号间隔使用更通用的界面，默认情况下使用Brent方法的一个版本。

• 第二个问题是能量的尺度是极端的，因此超出了求根方法的设计参数。一种解决方案是处理绝对和相对公差，以适应您的问题域和范围范围。另一种解决方案是将问题转换为具有更合理规模的域，以便错误控制方法/启发式在其设计范围和默认容差内工作。

sol = optimize.root_scalar(lambda s: eigen2(s*e_el), bracket=(100,200)) 
groundh = sol.root
groundh_ = sol.root*e_el

工作完美(可能使用布伦特方法作为括号方法的标准)，找到了138.023972 eV的基态能量，并得出了波形图

继续搜索，首先搜索符号更改，然后搜索根，在50*m*e_el到50*(m+1)*e_el的间隔内，查找下一个状态为

这篇关于量子谐振子/简谐振子模型的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！