python动力学建模与仿真_发动机悬置python仿真计算

代码可复制到https://hub.gke2.mybinder.org/user/lijil168-requirements-l6zexquh/tree运行

1、发动机悬置模态及解耦参考

2、发动机动力学激励计算参考

3、用数组和矩阵两种方式计算刚度矩阵，并对比结果，原文公式推导有点点错误。python用多维数组完成矩阵运算，很简洁而且可读性好。

4、思路：由发动机爆压计算倾覆力矩，建立状态空间动力学模型，进行模态及解耦率计算，最后进行仿真，得到发动机质心的状态，再换算到各个悬置点。

5、结果

分别为质量矩阵、刚度矩阵、自然频率、模态振型矩阵(每列对应振型向量)、解耦率(每列对应模态阶次，行从上到下对应x,y,z,theta_x,theta_y,theta_z 6个自由度)、阻尼矩阵。

5.1发动机激励计算结果

5.1悬置系统仿真结果

6、完整的仿真代码

​

# _*_ coding:UTF-8 _*_

import numpy as np

from numpy import sin

from numpy import cos

from numpy import arcsin as asin

from numpy import sqrt

from numpy import tan

from matplotlib.pyplot import*

from scipy import interpolate

def M_sys(m,Ixx,Ixy,Ixz,Iyy,Iyz,Izz):

return np.array([[m,0,0,0,0,0],\

[0,m,0,0,0,0],\

[0,0,m,0,0,0],\

[0,0,0,Ixx,-Ixy,-Ixz],\

[0,0,0,-Ixy,Iyy,-Iyz],\

[0,0,0,-Ixz,-Iyz,Izz]])

def K_element_i(k_pi,k_qi,k_ri,theta_pi,phi_pi,psi_pi,theta_qi,phi_qi,psi_qi,theta_ri,phi_ri,psi_ri):

#k_pi:第i个悬置的p向刚度；

#theta_pi：p轴与系统x轴的夹角；

#phi_pi：p轴与系统y轴的夹角；

#psi_pi:p轴与系统z轴的夹角；

K_xxi=k_pi*cos(theta_pi)**2+k_qi*cos(theta_qi)**2+k_ri*cos(theta_ri)**2

K_yyi=k_pi*cos(phi_pi)**2+k_qi*cos(phi_qi)**2+k_ri*cos(phi_ri)**2

K_zzi=k_pi*cos(psi_pi)**2+k_qi*cos(psi_qi)**2+k_ri*cos(psi_ri)**2

K_xyi=k_pi*cos(theta_pi)*cos(phi_pi)+k_qi*cos(theta_qi)*cos(phi_qi)+k_ri*cos(theta_ri)*cos(phi_ri)

K_xzi=k_pi*cos(theta_pi)*cos(psi_pi)+k_qi*cos(theta_qi)*cos(psi_qi)+k_ri*cos(theta_ri)*cos(psi_ri)

K_yzi=k_pi*cos(psi_pi)*cos(phi_pi)+k_qi*cos(psi_qi)*cos(phi_qi)+k_ri*cos(psi_ri)*cos(phi_ri)

return [K_xxi,K_xyi,K_xzi,K_yyi,K_yzi,K_zzi]

def K_sys(K_xxi,K_xyi,K_xzi,K_yyi,K_yzi,K_zzi,Ai,Bi,Ci):

Ai=np.array(Ai)

Bi=np.array(Bi)

Ci=np.array(Ci)

K_xx=sum(K_xxi)

K_xy=sum(K_xyi)

K_xz=sum(K_xzi)

K_yy=sum(K_yyi)

K_yz=sum(K_yzi)

K_zz=sum(K_zzi)

K_alpha_alpha=sum(K_yyi*Ci**2)+sum(K_zzi*Bi**2)-2*sum(K_yzi*Bi*Ci) ###

K_beta_beta=sum(K_xxi*Ci**2)+sum(K_zzi*Ai**2)-2*sum(K_xzi*Ai*Ci)

K_gamma_gamma=sum(K_xxi*Bi**2)+sum(K_yyi*Ai**2)-2*sum(K_xyi*Ai*Bi)   ####原公式有误！！！

K_alpha_beta=sum(K_xzi*Bi*Ci)+sum(K_yz*Ai*Ci)-sum(K_zzi*Ai*Bi)-sum(K_xyi*Ci**2)

K_beta_gamma=sum(K_xyi*Ai*Ci)+sum(K_xz*Ai*Bi)-sum(K_xxi*Ci*Bi)-sum(K_yzi*Ai**2) ###

K_alpha_gamma=sum(K_xyi*Ci*Bi)+sum(K_yz*Ai*Bi)-sum(K_yyi*Ci*Ai)-sum(K_xzi*Bi**2)

K_x_alpha=sum(K_xzi*Bi)-sum(K_xyi*Ci)

K_x_beta=sum(K_xxi*Ci)-sum(K_xzi*Ai)

K_x_gamma=sum(K_xyi*Ai)-sum(K_xxi*Bi)

K_y_alpha=sum(K_yzi*Bi)-sum(K_yyi*Ci)

K_y_beta=sum(K_xyi*Ci)-sum(K_yzi*Ai)

K_y_gamma=sum(K_yyi*Ai)-sum(K_xyi*Bi)

K_z_alpha=sum(K_zzi*Bi)-sum(K_yzi*Ci)

K_z_beta=sum(K_xzi*Ci)-sum(K_zzi*Ai)

K_z_gamma=sum(K_yzi*Ai)-sum(K_xzi*Bi)

K=np.array([[K_xx,K_xy,K_xz,K_x_alpha,K_x_beta,K_x_gamma],\

[K_xy,K_yy,K_yz,K_y_alpha,K_y_beta,K_y_gamma],\

[K_xz,K_yz,K_zz,K_z_alpha,K_z_beta,K_z_gamma],\

[K_x_alpha,K_y_alpha,K_z_alpha,K_alpha_alpha,K_