python发动机悬置解耦计算-按重心处整车坐标系解耦 一、参考二、python计算结果三、python程序四、virtual.lab motion模态计算验证五、Excite PU模态计算验证
一、参考
http://www.doc88.com/p-6834780959259.html
二、python计算结果
三、python程序 # _*_ coding:UTF-8 _*_import numpy as npfrom numpy import cosdef 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_sys_byMatrix(k_pi,k_qi,k_ri,theta_pi,phi_pi,psi_pi,theta_qi,phi_qi,psi_qi,theta_ri,phi_ri,psi_ri,Ai,Bi,Ci): #参数均采用数组格式 n=np.size(Ai) ki=np.zeros((n,3,3)) Ti=np.zeros_like(ki) BBi=np.zeros((n,3,6)) Ki=np.zeros((n,6,6)) #K=np.zeros((6,6)) for i in range(n): ki[i]=np.array([[k_pi[i],0,0], [0,k_qi[i],0], [0,0,k_ri[i]]]) Ti[i]=np.array([[cos(theta_pi[i]),cos(phi_pi[i]),cos(psi_pi[i])], [cos(theta_qi[i]),cos(phi_qi[i]),cos(psi_qi[i])], [cos(theta_ri[i]),cos(phi_ri[i]),cos(psi_ri[i])]]) BBi[i]=np.array([[1,0,0,0,Ci[i],-Bi[i]], [0,1,0,-Ci[i],0,Ai[i]], [0,0,1,Bi[i],-Ai[i],0]]) Ki[i]=BBi[i].T@Ti[i].T@ki[i]@Ti[i]@BBi[i] K=np.sum(Ki,axis=0) return Kdef Energy_DistributionJ(M,Value,Vector_Matrix): #3维能量分布：（阶次、Dof、Dof) n=M.shape[0] KE_klj=np.zeros((n,n,n)) for j in range(n): Value_j=Value[j] #0维数组：特征值 频率的平方 Vector_r=Vector_Matrix[:,j] #一维数组：行向量 Vector_c=Vector_r.reshape(n,1) #改成二维数组：列向量 KE_klj[j]=0.5*Value_j*M*Vector_c*Vector_r return KE_kljdef Energy_percent(KE_klj): n=KE_klj.shape[0] #KE_klj 3维能量分布：（阶次、Dof、Dof=页、行、列） #将每行Dof的能量合并缩维，得到（各Dof总能量占比,阶）： Energy_perDof=np.sum(KE_klj,2) #按行合并，第3维压缩掉，成2维数组：矩阵（阶次、各Dof能量） Energy_AllDof_r=np.sum(Energy_perDof,1) #按行合并，第2维压缩掉，成一维数组：行向量[1阶总能量、2阶总能量、...] Energy_AllDof_c=Energy_AllDof_r.reshape(n,1) #改成二维数组：列向量[[1阶总能量],[2阶总能量]、...] Energy_percent=100.0*Energy_perDof/Energy_AllDof_c #2维数组：矩阵（阶次、各Dof能量占比） Energy_percent=Energy_percent.T ########二维数组转置，得到矩阵（各Dof总能量占比,阶） return Energy_percent###########input###########################################################################[m,Ixx,Iyy,Izz,Ixy,Ixz,Iyz]=[1224.4,79.385,523.73,491.332,-5.866,-76.403,-4.698]PG=[0.266249,-0.011964,0.186296]ABC=[[-0.525,-0.413,0.110],[-0.525,0.413,0.110],[0.675,-0.380,0.151],[0.675,0.380,0.151]]ABC=np.array(ABC)-np.array(PG)Ai,Bi,Ci=ABC[:,0],ABC[:,1],ABC[:,2]theta_pi=[0*np.pi/180,0*np.pi/180,0*np.pi/180,0*np.pi/180] #3theta_qi=[90*np.pi/180,90*np.pi/180,90*np.pi/180,90*np.pi/180]theta_ri=[90*np.pi/180,90*np.pi/180,90*np.pi/180,90*np.pi/180] #93phi_pi=[90*np.pi/180,90*np.pi/180,90*np.pi/180,90*np.pi/180]phi_qi=[0,0,0,0]phi_ri=[90*np.pi/180,90*np.pi/180,90*np.pi/180,90*np.pi/180]psi_pi=[90*np.pi/180,90*np.pi/180,90*np.pi/180,90*np.pi/180] #87psi_qi=[90*np.pi/180,90*np.pi/180,90*np.pi/180,90*np.pi/180]psi_ri=[0*np.pi/180,0*np.pi/180,0*np.pi/180,0*np.pi/180] #3k_pi=[1386000,1386000,2745000,2745000]k_qi=[354000,354000,367500,367500]k_ri=[1068000,1068000,2025000,2025000]###########end input###########################################################################M_sys=M_sys(m,Ixx,Ixy,Ixz,Iyy,Iyz,Izz)#[K_xxi,K_xyi,K_xzi,K_yyi,K_yzi,K_zzi]=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=K_sys(K_xxi,K_xyi,K_xzi,K_yyi,K_yzi,K_zzi,Ai,Bi,Ci)#print(K*(abs(K)>0.001))K_sys=K_sys_byMatrix(k_pi,k_qi,k_ri,theta_pi,phi_pi,psi_pi,theta_qi,phi_qi,psi_qi,theta_ri,phi_ri,psi_ri,Ai,Bi,Ci)Value, Vector_Matrix = np.linalg.eig(np.linalg.inv(M_sys)@K_sys)#np.set_printoptions(formatter={'float': '{: 0.3f}'.format})np.set_printoptions(precision=3, suppress=True)print(M_sys)print(K_sys)print(np.sqrt(Value)/2/np.pi)print(Vector_Matrix)#print(Vector_Matrix.T@M_sys@Vector_Matrix)#print(Vector_Matrix.T@K_sys@Vector_Matrix)KE_klj=Energy_DistributionJ(M_sys,Value,Vector_Matrix)print(Energy_percent(KE_klj)) 四、virtual.lab motion模态计算验证
五、Excite PU模态计算验证
当地板水平放置，结果见上面，当地板逆时针转45度，模态结果见下图，自然频率不变。
virtual.lab motion （使用BUSH力单元）和PU（使用NONL）模态计算的结果完全一致，和理论计算接近。但是如果PU使用emo1，在地板偏转一个角度，自然频率会出现较大的变化，原因不明！！
见下图：