概述
我想找到微分方程的平衡点,并检查平衡点是否稳定.
这是一个最小的工作示例
import numpy as np
from scipy.optimize import fsolve
dim = 2
A = np.random.uniform(size = (dim,dim))
sol,infodict,ier,mesg = fsolve(lambda x: 1-np.dot(A,x),np.ones(dim),full_output = True)
为了找出解sol是否稳定,雅可比行列的所有特征值必须具有负实部.然而,Jacobian没有保存在infodict中,而是QR分解被保存在infodict中.
如何从fsolve的QR分解中获得Jacoian？
我能做的就像是
eigenvalues = np.linalg.eigvals(infodict["fjac"])*infodict["r"][ind]
ind是r的对角线条目,但我怀疑这是最好的方法.
结论是通过乘法重建雅可比行列是这里的方法.你在做什么(仅仅找到Q因子的特征值？)是不正确的.首先,使用np.triu_indices从给定的平面形式恢复R矩阵;然后将Q乘以R;然后找到特征值.
r = np.zeros((dim,dim))
r[np.triu_indices(dim)] = infodict["r"]
eigenvalues = np.linalg.eigvals(infodict["fjac"].dot(r))
总结
以上是编程之家为你收集整理的使用fsolve检查微分方程的稳定性全部内容，希望文章能够帮你解决使用fsolve检查微分方程的稳定性所遇到的程序开发问题。
如果您也喜欢它,动动您的小指点个赞吧
602392714
清零编程群