0

我正在尝试计算一个矩阵的特征值F，它包含 9 个变量，它们是三个向量的笛卡尔坐标。的执行时间F.eigenvals()至少需要 15 分钟，而且我没有足够的耐心等待程序完成更长时间。我的怀疑是我做错了可怕的事情，这导致了很长的计算。也许是 numpy 和 sympy 对象的混合导致了这个问题。

M我有三个水原子和它们各自的质量，我使用函数计算一个所谓的质量矩阵Mmat。

import sympy as sy
import numpy as np
from scipy.constants import u as u

def Mmat(*args):

    N = len(args)
    args = np.array([[arg*u]*3 for arg in args]).reshape(N*3)
    m = args**(-1/2)
    I = np.identity(3*N)

    return m*I

mO, mH1, mH2 = 16, 1, 1
M = Mmat(mH1,mH2,mO)

该质量矩阵需要根据公式F = M*H*M与 Hessian 矩阵相乘，H以生成 (9,9) 矩阵F，从中提取特征值。Hessian 矩阵本身由V三个原子 O、H1 和 H2 的笛卡尔坐标的势的偏导数定义。

矩阵F本身的计算只需要一秒钟，但应用该.eigenvals()方法需要很长时间。这是计算的函数H：

def dist(v1,v2):

    return ((v2-v1)**2).sum()**(-1/2)

def Hessian():

    kOH, kHH, dr1, dr2, dr3 = sy.symbols('kOH kHH dr1 dr2 dr3', real=True, positive=True)

    V = 1/2*kOH*(dr1)**2 +1/2*kOH*(dr2)**2 +1/2*kHH*(dr3)**2

    x1,x2,x3,x4,x5,x6,x7,x8,x9 = sy.symbols('x1 x2 x3 x4 x5 x6 x7 x8 x9', real=True)
    dOH10, dOH20, dH1H10 = sy.symbols('dOH10 dOH20 dH1H10', real=True, positive=True)

    vO = np.array([x1,x2,x3])
    vH1 = np.array([x4,x5,x6])
    vH2 = np.array([x7,x8,x9])

    r1 = dist(vO,vH1) - dOH10
    r2 = dist(vO,vH2) - dOH20
    r3 = dist(vH1,vH2) - dH1H10

    V = V.subs(dr1,r1).subs(dr2,r2).subs(dr3,r3)

    H = sy.hessian(V,[x1,x2,x3,x4,x5,x6,x7,x8,x9])

    return H

F = M*Hessian()*M

真正的问题只出现在最后：

print(F.eigenvals())

我是否不恰当地混合了对象类？质量矩阵（6.022e26）中有大量数字是事实吗？这种大小的矩阵运行如此缓慢是否正常？

先感谢您。

4

0 回答 0