计算物理专题:有限差分法解决本征值问题

news2025/1/17 0:07:46

计算物理专题:有限差分法解决本征值问题

定态薛定谔方程差分形式

  • 一维定态薛定谔方程

-\alpha\frac{d^2}{dx^2}\psi(x)+V(x)\psi(x)=E\psi(x);\alpha=\frac{\hbar^2}{2\mu}

\left\{\begin{matrix} \psi(kh)=\psi_k\\ R=\alpha/h^2\\ a_k=2R+V(kh) \end{matrix}\right.

\begin{bmatrix} a_1& -R & & & \\ -R&a_2 &-R & & \\ & -R& a_3 &-R & \\ & & & \ddots& \\ & & & -R & a_n \end{bmatrix}\begin{bmatrix} \psi_1\\ \psi_2\\ \psi_3\\ \cdots\\ \psi_n \end{bmatrix}=E\begin{bmatrix} \psi_1\\ \psi_2\\ \psi_3\\ \cdots\\ \psi_n \end{bmatrix}

谐振子

V(x)=\frac{1}{2}\hbar\omega^2x^2

\hbar=1;\omega=1;\mu=1

x_0=0;x_1=20;h=0.01;

解法代码

import numpy as np
def householder(symmetric_matrix):
    M = symmetric_matrix
    assert np.allclose(M,M.T),"matrix is not symmetric"
    N = len(M)
    for i in range(1,N-1):
        r = 0
        for j in range(i,N):
            r += M[i-1][j]**2
        r = r**0.5
        if r * M[i-1][i] > 0:
            r *= -1
        Ki = -1/(r**2 - r*M[i-1][i])
        Pi = np.zeros((N,1))
        Pi[i][0] = M[i-1][i] - r
        for j in range(i+2,N+1):
            Pi[j-1][0] = M[i-1][j-1]
        Fi = np.dot(Pi,Pi.T)*Ki
        for j in range(1,N+1):
            Fi[j-1][j-1] += 1        
        M = np.dot(np.dot(Fi,M),Fi)
    return M

def qr_decomposition(matrix):
    m,n = matrix.shape
    Q = np.zeros((m,n))
    R = np.zeros((n,n))
    for j in range(n):
        v = matrix[:,j]
        for i in range(j):
            R[i,j] = np.dot(Q[:,i],matrix[:,j])
            v = v - R[i,j] * Q[:,i]
        R[j,j] = np.linalg.norm(v)
        Q[:,j] = v/R[j,j]
    return Q,R

def QLmethod(matrix,tol=1e-6,maxiter=100):
    assert np.allclose(matrix,matrix.T),"matrix is not symmetric"
    matrix = householder(matrix)
    n = matrix.shape[0]
    eigenvalues = np.zeros(n)
    iterations = 0
    while np.max(np.abs(matrix.diagonal(offset=1))) > tol and iterations < maxiter:
        q,r = qr_decomposition(matrix)
        matrix = np.dot(r,q)
        iterations += 1
    eigenvalues = matrix.diagonal()
    return eigenvalues


hbar = 1
omega = 1
mu = 1
x0 = 0
x1 = 10
h = 0.1

V = lambda x:1/2*hbar*omega**2*x**2
R = hbar**2/2/mu/h**2
N = int((x1-x0)/h)
List = [2*R + V(k*h) for k in range(1,1+N)]

symmetric_matrix = np.diag(List)
for i in range(N-1):
    symmetric_matrix[i][i+1] = -R
    symmetric_matrix[i+1][i] = -R

#Three Symmetric matrix
eigenvalue = QLmethod(symmetric_matrix)

运行结果

array([199.05965221, 197.85154038, 196.69647447, 195.59251125,
       194.53779897, 193.53057225, 192.56914726, 191.65191722,
       190.7773482 , 189.94397515, 189.15039823, 188.39527932,
       187.6773388 , 186.99535243, 186.34814856, 185.73460537,
       185.15364836, 184.60424796, 184.08541726, 183.59620995,
       183.13571823, 182.70307101, 182.29743212, 181.91799858,
       181.56399908, 181.23469247, 180.92936632, 180.64733559,
       180.38794141, 180.15054982, 179.93455068, 179.73935657,
       179.56440178, 179.40914085, 179.27304046, 179.15548747,
       179.0549045 , 178.96231954, 178.82813537, 178.45299393,
       177.32375035, 174.73750414, 170.45930442, 165.03198466,
       159.18519843, 153.34588073, 147.67451019, 142.21724325,
       136.98169564, 131.96220481, 127.1481637 , 122.52724056,
       118.08680316, 113.81459325, 109.69905122, 105.72946646,
       101.89603735,  98.18988381,  94.60303514,  91.12840579,
        87.7597663 ,  84.49171355,  81.31964293,  78.2397245 ,
        75.24888443,  72.34479386,  69.52586706,  66.791272  ,
        64.14095638,  61.57569128,  59.0971279 ,  56.70783771,
        54.41122812,  52.21100812,  50.10939318,  48.1026652 ,
        46.17390347,  44.28931575,  42.41015044,  40.51421357,
        38.60050745,  36.67629929,  34.74758658,  32.81748476,
        30.88705093,  28.95607702,  27.02366776,  25.08868837,
        23.15008972,  21.20709216,  19.25922869,  17.30628391,
        15.34819203,  13.38495142,  11.41657896,   9.44309464,
         7.46451849,   5.4808703 ,   3.49216962,   1.49843574])

存疑

  • 谐振子势确实是等距的,但是为什么计算出的结果与E_n=(n+\frac{1}{2})\hbar\omega 相差了一个1,这是为什么呢?没想明白。

向Fortran的翻译

有限差分法解决本征值问题的特点

  • 引用有限差分法解决本征值问题时,步长的大小和边值条件会极大得影响计算的精度,有限差分法一般可以提供很好的近似解

量子隧穿效应的调研

\subsubsection{量子隧穿与核聚变}
    恒星原子核之间存在很高的排斥静电力。对于核聚变,需要大量的能量使原子核彼此足够接近以克服库仑力:
    原子核的平均动能可以被表示为
    \begin{equation}
        E_{kin}=\frac{3}{2}kT
    \end{equation}
    需要被克服的能量$E_{cb}$为:
    \begin{equation}
        E_{cd}=\frac{Z_1Z_2e^2}{r_c}
    \end{equation}
    由于恒星内部的原子核在静止燃烧过程中无法克服库伦势垒,因此除了热运动以外,必须要存在其他因素才能使得热核反应成为可能。尽管恒星内部的热运动低到无法直接引起热核反应,但它能在库仑势垒后面产生显著的概率密度使得量子隧穿能以一定的概率发生。这一事件所需要的能量远远低于克服库仑势垒所需的动能。对于运动速度为$\nu$的核子发生量子隧穿的概率$P_t(\nu)$正比于:
    \begin{equation}
        \exp{\frac{-4\pi ^2Z_1 Z_2 e^2}{h}\frac{1}{\nu}}
    \end{equation}

\subsubsection{分子生物学中量子隧穿}
    比氢原子重的化学元素合成与量子隧穿现象密切相关。对氘代细菌孢子的实验表明,与水基培养基中的自发突变率相比,在含氘氧化物的培养基中生长的细胞中的自发点突变率较低。这意味着突变偏差对所涉及的粒子质量的依赖性。而量子隧穿速率强烈依赖于粒子质量。$Löwdin DNA  $突变模型描述了互补碱基对之间氢键所涉及的质子的量子性质。$Löwdin$的模型基于稳定的沃森-克里克碱基配对和稳定的遗传信息需要质子不改变其在碱基对中的位置这一前提,将质子在碱基对中的位置用所涉及的核酸碱基分子的规范形式来确定。$Löwdin$模型指出,通过质子对转移势垒的量子隧穿,当质子在DNA复制后立即处于量子力学非稳态时,质子在核酸碱基内的位置可能会以很小的概率因量子隧穿发生变化。\par 
    Lödwin模型描述了DNA中质子的量子隧穿导致的普遍点突变偏差。更多的研究表明量子隧穿不仅能起到对DNA的破坏,而且还涉及到DNA修复问题。紫外线照射可以通过诱导DNA链中凸起的形成来损伤DNA。这些凸起是由于DNA链上相邻的嘧啶在紫外线照射下发生二聚而引起的。黄蛋白光解酶能够通过电子转移分裂共价连接的嘧啶来修复这种变形的DNA,从而使它们恢复到正常的单体形式,不再形成二聚体。长距离的超交换介导隧穿使得蛋白质中电子的量子隧穿可以在高达3纳米的距离上发生,这种势垒宽度远大于真空中的可能宽度。蛋白质通过超能变化介导的机制提供的低电子态提高了电子转移速率,即使在如此长的距离上它们仍然相对较高。蛋白质对损伤DNA的修复基本上基于长程电子隧穿。\par 
    

\subsubsection{量子隧穿与语音识别$^{[3]}$}
    语音信号通常服从两种分布:较短语音信号服从Gauss分布,较长语音信号服从Laplace分布。对于同一说话人,由于自身身体条件, 对于同样的语言内容,发音频率总在误差允许的范围内相似。对于不同的说话人,一方面发出的声音信号存在差异,这种差异表现为频率特征上的差异,可由处于不同稳定的量子态来描述。\par 
    另一方面,对不同的语音信号分帧处理后,由于没帧的时间很短,基本上服从Gauss分布。因此,每一分帧语音信号可以被视作一个包含一组频率特征的量子波函数。不同的频率,通过相同的势垒,其隧穿系数不同,隧穿后的频率不同。\par 
    因此,通过设置一组势垒,让每一个势垒有唯一的一个频率透射即可构成一组特征向量。这些特征都是非负的,再构成一个随机向量,可以用高斯向量表征。根据向量中元素拟合,降维后成的二维概率密度函数,通过最大似然估计就能实现对说话人的识别。\par 

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/678496.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

chatgpt赋能python:用Python分析电影评分数据

用Python分析电影评分数据 Python是一种流行的数据分析和可视化工具&#xff0c;它可以让我们更深入地了解电影的评分数据。在本文中&#xff0c;我们将使用Python来分析一些电影评分数据&#xff0c;并试图找出一些有趣的模式和趋势。 数据来源 我们将使用公共数据集IMDb电…

第4章 网络层

1‌、下列关于路由算法描述错误的是&#xff08; &#xff09; A. 链路状态算法是一种全局路由算法&#xff0c;每个路由器需要维护全局状态信息B. OSPF 是一种域内路由协议&#xff0c;核心是基于 Dijkstra 最低费用路径算法C. RIP 是一种域内路由算法&#xff0c;核心是基…

采用SqlSugar的DBFirst相关功能创建数据库表对应的实体类

.NET Core官方教程中推荐使用的EF Core数据库ORM框架虽然能用&#xff0c;但是用起来并不是太方便&#xff08;或者是不习惯&#xff0c;之前用的最多的还是linq&#xff09;。之前下载的开源博客项目中使用的SqlSugar&#xff0c;后者是由果糖大数据科技团队维护和更新 &#…

基于WebAssembly构建Web端音视频通话引擎

Web技术在发展&#xff0c;音视频通话需求在演进&#xff0c;怎么去实现新的Web技术点在实际应用中的值&#xff0c;以及给我们带来更大的收益是需要我们去探索和实践的。LiveVideoStackCon 2022北京站邀请到田建华为我们从实践中来介绍WebAssembly、WebCodecs、WebTransport等…

【裸机开发】IRQ 中断服务函数(一) —— 汇编初始化

IRQ 和前面的Reset 函数不大一样&#xff0c;当一个IRQ中断产生时&#xff0c;我们也不知道这个IRQ中断来自哪个外设&#xff0c;因此&#xff0c;需要先获取到中断ID&#xff0c;随后才会跳转到真正的中断服务函数执行处理逻辑。 整个 IRQ 中断处理可以看做是包含了两个部分&…

MySQL 自增主键一定是连续的吗?

众所周知&#xff0c;自增主键可以让聚集索引尽量地保持递增顺序插入&#xff0c;避免了随机查询&#xff0c;从而提高了查询效率 但实际上&#xff0c;MySQL 的自增主键并不能保证一定是连续递增的。 下面举个例子来看下&#xff0c;如下所示创建一张表&#xff1a; 自增值保…

ORCA优化器浅析——GP数据库调用优化器流程

首先我们需要看CGPOptimizer类(src/include/gpopt/CGPOptimizer.h)为Greenplum数据库提供ORCA优化器export出来的函数的封装。Greenplum数据库主流程调用extern "C"中提供的函数&#xff0c;比如初始化ORCA优化器的函数InitGPOPT&#xff0c;优化查询树的函数GPOPTOp…

springboot+jsp农产品商城宣传网站设计与实现oo6e3

在该在线助农系统设计与实现中&#xff0c;idea能给用户提供更多的方便&#xff0c;其特点一是方便学习&#xff0c;方便快捷&#xff1b;二是有非常大的信息储存量&#xff0c;主要功能是用在对数据库中查询和编程。其功能有比较灵活的数据应用&#xff0c;只需利用小部分代码…

【Leetcode60天带刷】day30回溯算法——332.重新安排行程 , 51. N皇后 ,37. 解数独

​ 题目&#xff1a; 332. 重新安排行程 给你一份航线列表 tickets &#xff0c;其中 tickets[i] [fromi, toi] 表示飞机出发和降落的机场地点。请你对该行程进行重新规划排序。 所有这些机票都属于一个从 JFK&#xff08;肯尼迪国际机场&#xff09;出发的先生&#xff0c;…

【从零开始学习JAVA | 第十四篇】继承

目录 前言&#xff1a; 引入&#xff1a; 继承&#xff1a; 小拓展&#xff1a; 优点&#xff1a; 成员方法的继承问题&#xff1a; 总结&#xff1a; 前言&#xff1a; 继承是面向对象三大特性之一&#xff0c;它是在封装之后我们讲解的一个重要的性质&#xff0c;继承…

在github上创建个人主页的方法【2023更新版】

01-进入github的网站&#xff0c;链接 https://github.com/ &#xff0c;然后注册&#xff0c;登陆&#xff0c;注意登陆时设置的用户名(username)就是将来你个人主页的三级域名&#xff0c;所以这里一定要慎重填写username。如下图所示&#xff1a; 02-注册完成后进入个人主…

2024考研408-计算机组成原理第四章-指令系统

文章目录 前言一、指令系统现代计算机的结构1.1、指令格式1.1.1、指令的定义1.1.2、指令格式1.1.3、指令—按照地址码数量分类①零地址指令②一地址指令&#xff08;1个操作数、2个操作数情况&#xff09;③二地址指令④三地址指令⑤四地址指令 1.1.4、指令-按照指令长度分类1.…

【计算机组成原理】Yy-z02模型机的硬布线控制器设计

目录 一、Yy-z02模型机的系统结构 二、Yy-z02模型机的数据通路 三、Yy-z02模型机的指令执行 四、Yy-z02模型机的硬布线控制器 一、Yy-z02模型机的系统结构 指令系统的实现 <--- 构造它的硬件系统 硬件系统构造过程&#xff1a; 分析指令格式和各指令的功能确定部件连…

金蝶软件遭遇.locked勒索病毒攻击:如何保护与解救您的数据?

引言&#xff1a; 近期&#xff0c;部分运行金蝶云星空软件的服务器遭受了一场勒索病毒的网络安全攻击&#xff0c;其重要数据遭到了.locked勒索病毒的加密。作为一个知名的企业级ERP软件及财务软件&#xff0c;金蝶软件的数据安全事关客户和企业的利益。91数据恢复在本文将深…

【王道·操作系统】第四章 文件管理【未完】

一、初识文件管理 文件&#xff1a;一组有意义的信息/数据集合文件属性&#xff1a; 文件名&#xff1a;创建文件的用户决定&#xff0c;主要是为了方便用户找到文件&#xff0c;同一目录下不允许有重名文件标识符&#xff1a;一个系统内的各文件标识符唯一&#xff0c;对用户来…

老大给了个新需求:如何将汉字转换成拼音字母?1行Python代码搞定!

大家好&#xff0c;这里是程序员晚枫&#xff0c;小红薯也叫这个名。 之前的视频给大家分享了&#xff1a;中文编程&#xff0c;一行代码实现。 今天给大家分享一下&#xff0c;如何通过1行Python代码&#xff0c;实现汉语转拼音 1、先上代码 实现汉语转拼音效果的第三方库…

逻辑回归(Logistics Regression)

1.逻辑回归&#xff08;Logistics Regression&#xff09; 逻辑回归用于解决二分类问题 1.1 Sigmoid函数 sigmoid函数在神经网络中如何起作用&#xff1f;详见本人笔记&#xff1a;机器学习和AI底层逻辑 复杂非线性分类->多个线段->每个线段是叠加而来的->sigmoid函…

计算机视觉 + Self-Supervised Learning 五种算法原理解析

计算机视觉领域下自监督学习方法原理 导语为什么在计算机视觉领域中进行自我监督学习&#xff1f; 自监督学习方法Generative methodsBEiT 架构 Predictive methodsContrastive methodsBootstraping methodsSimply Extra Regularization methods 导语 自监督学习是一种机器学习…

Web服务器群集:Nginx网页及安全优化

目录 一、理论 1.Nginx网页优化 2.Nginx安全优化 3.Nginx日志分割 二、实验 1.网页压缩 2.网页缓存 3.连接超时设置 4.并发设置 5.隐藏版本信息 6.脚本实现每月1号进行日志分割 7.防盗链 三、总结 一、理论 1.Nginx网页优化 &#xff08;1&#xff09;概述 在企…

神仙打架的618,谁才是真正的大赢家?

618大促已经缓缓落下帷幕&#xff0c;各大平台和品牌方准时准点晒出成绩单。在一串又一串红彤彤的战报中&#xff0c;家电品牌你追我赶的激烈战况一如以往。 我们从中也得以窥见新消费时代下中国家电行业的未来&#xff0c;尤其是在消费者纷纷捂紧钱袋子的今年&#xff0c;红色…