计算物理专题：实对称矩阵特征值的求解问题 |
您所在的位置：网站首页 › Python矩阵特征值 › 计算物理专题：实对称矩阵特征值的求解问题 |
计算物理专题：实对称矩阵特征值的求解问题
|
我们需要求出一个实对称阵的全部特征值与特征向量实对称矩阵特征值的求解是计算量子力学的基础之一（算符是厄米的）同时，实对称矩阵的求解也是使用有限差分法解本征方程的重要方法之一
雅克比方法
我知道这个推导过程非常的复杂，所以我不想推导在这里，确实要明白的，可以私聊
雅克比方法是求解实对称矩阵的全部本征值和本征矢量的最常用方法雅克比方法虽然能求出全部的本征值，但是在阶数较高的矩阵上计算复杂度很高
import numpy as np
import math
def sign(x):
if x>0:
return 1.0
elif x= max_iter:
print("end!!")
break
Lambda = abs(A[j,j]-A[i,i])
Mu = 2.0 * A[i,j] * sign(A[j,j]-A[i,i])
if Lambda**2+Mu**2 != 0:
Omega = Lambda/(Lambda**2+Mu**2)**0.5
c = ((1+Omega)/2)**0.5
s = Mu*Omega/(Lambda*(2*(1+Omega))**0.5)
else:
c = 2**0.5/2
s = 2**0.5/2
J = np.eye(n)
J[i, i] = c
J[j, j] = c
J[i, j] = s
J[j, i] = -s
A = np.dot(np.dot(J.T, A), J)
V = np.dot(V, J)
iterations += 1
eigenvalues = np.diag(A)
eigenvectors = V.T
return eigenvalues, eigenvectors
A = np.array([[1,2,1],
[2,4,0],
[1,0,3]])
eigenvalues, eigenvectors = jacobi_eigenvalue(A)
print("本征值:")
print(eigenvalues)
print("本征向量:")
print(eigenvectors)
QL方法
首先利用householder 变换（不改变矩阵的特征值）将实对称矩阵变换成一个对称三对角矩阵然后利用QL方法求解这个对称三对角矩阵的特征值
QL方法的计算速度会更快一些
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
symmetric_matrix = np.array([[1,-7,-2,1],
[-7,40,-13,26],
[-2,-13,3,-12],
[1,26,-12,71]],dtype="float64")
#Three Symmetric matrix
eigenvalue = QLmethod(symmetric_matrix)
经典电动力学是物理学专业基础课的重要分支，它研究电荷和电磁场之间的相互作用规律。电动力学理论不仅是现代物理学的基础，而且在电子学、通信、能源、材料科学、生物医学等领域都有广泛的应用。 本次摘要笔记旨在简述，评注，补充经典电动力学的基本概念、基本定律和应用。笔记内容包括电荷、电场、磁场、电流、电磁波等基本概念，以及库仑定律、电场的高斯定理、电场的环路定理、法拉第电磁感应定律和安培定律等基本定律。此外，笔记还简要得总结了绝大多数课堂报告的主要内容。课堂报告的主要内容已经上传到了我的博客。 除了阅读郭硕鸿老师的教材，我还阅读了John David Jackson的经典教材。最终我阅读到了chapter 12 Proca Lagrangian；Photon Mass Effects。这本书补充了许多的内容，难度也非常的大，总体上还是非常有趣的。我读得很快，很多地方只是囫囵吞枣，了解个大概。 大学二年级下学期的学习并不是非常的顺利。在大二上学期学习完了数理方法和抽象代数之后，已经可以从许许多多的方面认识新的知识了，这是一种非常奇妙的体验。 这学期修读的压力还是很大的，很多时候我都在怀疑自己是否能够很好的在物理学上做下去。我并不像是一个有着很好的数学能力和记忆力的孩子或者真的能理解每一门课程的核心知识的学生。但是我想也许以后还有很多很多的选择和机会吧。 这学期，我也了解了计算电磁学。借助计算物理的机会，也练习了自己。 我认识一个非常厉害外院的同学，他一直想修读物理专业，两年来几乎选修了物理系每一门专业课并都取得了非常优秀的成绩。但是因为种种原因，他并不能如愿转到物理系来。生活就是这么神奇吧。我记得他说过的一句话：“无论在不在物理界待着，都要运用理性的思维严密的逻辑去解决生活的问题，然后积极的迎接新的一天。”如果一定要说世界上有什么珍贵的东西的话，那么我想也无非自由漫步的思想，探穷无尽的自然与值得依靠的亲友吧。 自适应光学的发展展望激光通信可能解决不适宜铺设光纤的场景的通信问题。大气中温差差，风速都可能导致频繁的湍流活动，从而诱发波前畸变，不利于光信号与单模光纤的耦合。配合上自适应光学，激光通信的准确性，通信长度被大大增强了。 从军事角度上看，城市中对传统的电磁信号干扰大，干扰多，传统的通信方法不能适应现代战争尤其是防空作战的需要，因而激光通信对城市防空有着重要的意义。自适应光学的发展提高了激光通信的效率与准确性。 作为基础学科，天文观测能力的提升意味着更多的极端物理环境可以被观测到，为验证物理学的许多理论提供了实验的验证。 下表是从万方主题分析报告检索关键词----自适应光学 + 激光通信 得到的数据：
|
【本文地址】
今日新闻 |
推荐新闻 |