单自由度系统的随机振动计算与abaqus仿真验证

单自由度系统的随机振动计算与abaqu仿真验证

Abaqus使用文档网页版

注意:Abaqus的很多关键字Abaqus CAE都无法识别,可以用inp文件进行处理,比如随机噪声激励力的加载,受迫响应激励表的输入等等.

Abaqus关键字参考key

Abaqus阻尼设置参考cdsn

1、参考单自由度振动计算

2、习题

假设一自由度系统受到的激励力是均值为0、自谱为S0=1的理想白噪声平稳过程。求系统响应的均方值与自谱（功率谱密度函数，表示随机过程的均方值在频率内的分布密度）。
设m=1,k=(2*pi)**2,阻尼系数zeta=0.16,即构造固有频率为1Hz,快速衰减的系统。
激励力随机时，和确定时类似，都需要利用系统的脉冲响应函数h(t)或者复频响应函数H(w)。
系统的脉冲响应函数h(t)见下图：

系统频响函数见下图：

理想白噪声激励力的自谱见下图：

位移响应的自谱见下图：

3、python计算代码

import matplotlib.pyplot as plt
import numpy as np

from numpy import pi,sqrt,cos,sin,arctan,exp

plt.rcParams['font.sans-serif']=['SimHei'] #用来正常显示中文标签
plt.rcParams['axes.unicode_minus']=False #用来正常显示负号

##########begin:计算单自由度振动###############
k=(2*pi)**2
m=1
wn=sqrt(k/m)
s_=0.02 #0.02 0.16
wd=sqrt(1-s_**2)*wn
c=2*s_*wn*m
delta=s_*wn

#系统激励力为均值为0自谱为S0的理想白噪声平稳过程：
S0=1
#系统的脉冲响应函数：
h_t=lambda t:1/(m*wd)*exp(-delta*t)*sin(wd*t)
#系统的频响函数：
H_w=lambda w:1/sqrt((k-m*w**2)**2+c**2*w**2)
#响应的均方值：
phi2_x=S0/(2*c*k)
#响应的自谱：
Sxx_w=lambda w:S0/((k-m*w**2)**2+c**2*w**2)

ti=np.array(np.linspace(0,4,1000,False))#0:0.002:2;
h_ti=h_t(ti)
wi=np.array(np.linspace(0,100,1000,False))#0:0.002:2;
H_wi=H_w(wi)
Sxx_wi=Sxx_w(wi) #双边功率谱的一边
Gxx_wi=2*Sxx_wi #单边功率谱

plt.figure
fig_mpl, ax =plt.subplots(4,1,figsize=(5,3), facecolor='white')
ax[0].plot(ti,h_ti,label='系统的脉冲响应函数')
ax[0].set_title("h(t)")
ax[0].set_xlabel("t")
ax[0].legend()
ax[1].plot(wi,H_wi,label='系统频响函数')
ax[1].set_title("H(w)")
ax[1].set_xlabel("w")
ax[1].legend()
ax[2].plot(wi,np.ones(np.size(wi))*S0*2,label='输入力白噪声自谱（单边）')
ax[2].set_title("Gff(w)")
ax[2].set_xlabel("w")
ax[2].legend() #指定legend
ax[3].plot(wi,Gxx_wi,label='输出位移自谱（单边）；均方值:%f'%phi2_x)
ax[3].set_title("Gxx(w)")
ax[3].set_xlabel("w")
plt.legend() #plt.legend只标注最后一张图
plt.subplots_adjust(left=None, bottom=None, right=None,
top=None, wspace=1, hspace=1)
plt.show()

4、Abaqus仿真验证

用Abaqus建立一自由度振动系统,两个part 点,用弹簧相连.受迫响应如下图.固有频率1Hz


inp文件:

*Heading
** Job name: Job-1 Model name: Model-1
** Generated by: Abaqus/CAE 2016
*Preprint, echo=NO, model=NO, history=NO, contact=NO
**
** PARTS
**
*Part, name=Part-1
*Node
      1,           0.,         100.,           0.
*Nset, nset=Part-1-RefPt_, internal
1, 
*End Part
**  
*Part, name=Part-2
*Node
      1,           0.,           0.,           0.
*Nset, nset=Part-2-RefPt_, internal
1, 
*End Part
**  
**
** ASSEMBLY
**
*Assembly, name=Assembly
**  
*Instance, name=Part-1-1, part=Part-1
*End Instance
**  
*Instance, name=Part-2-1, part=Part-2
*End Instance
**  
*Nset, nset=Set-1, instance=Part-1-1
 1,
*Nset, nset=_PickedSet29, internal, instance=Part-2-1
 1,
*Nset, nset=_PickedSet32, internal, instance=Part-1-1
 1,
*Element, type=MASS, elset=_PickedSet31_Inertia-1_
1, Part-1-1.1
*Mass, elset=_PickedSet31_Inertia-1_
1., 
*Spring, elset=Springs/Dashpots-1-spring

39.4384
*Element, type=SpringA, elset=Springs/Dashpots-1-spring
2, Part-1-1.1, Part-2-1.1
*End Assembly
** ----------------------------------------------------------------
** 
** Force Excitation
** 
*AMPLITUDE, NAME=Powerunit-3376-2_R, DEFINITION=TABULAR, VALUE=RELATIVE,
INPUT=Powerunit-3376-2_R.cload
*AMPLITUDE, NAME=Powerunit-3376-2_I, DEFINITION=TABULAR, VALUE=RELATIVE,
INPUT=Powerunit-3376-2_I.cload
**
**
*Psd-definition, name=Amp-1, g=9800.
1., 0., 0.1
1., 0., 100.
** 
** BOUNDARY CONDITIONS
** 
** Name: BC-1 Type: Symmetry/Antisymmetry/Encastre
*Boundary
_PickedSet29, ENCASTRE
** Name: BC-2 Type: Symmetry/Antisymmetry/Encastre
*Boundary
**_PickedSet32, YASYMM
_PickedSet32,1,1,0
_PickedSet32,3,6,0
** ----------------------------------------------------------------
** 
** STEP: Step-1-modal
** 
*Step, name=Step-1-modal, nlgeom=NO, perturbation
*Frequency, eigensolver=Lanczos, sim=NO, acoustic coupling=on, normalization=displacement
, , 1000., , , 
** 
** OUTPUT REQUESTS
** 
*Restart, write, frequency=0
** 
** FIELD OUTPUT: F-Output-1
** 
*Output, field
*Node Output, nset=Set-1
UT, 
*End Step
** ----------------------------------------------------------------
** 
** STEP: Step-2-Excite FEA Dynamic Reponse Analysis
** 
*STEP
*STEADY STATE DYNAMICS, FREQUENCY SCALE=LINEAR, INTERVAL=RANGE
0.00000,29.9147726,118,1
** 
** 
** Modal Damping
** 
*MODAL DAMPING, RAYLEIGH
1,36,41.8879020479,1.48544613552e-005
** 
*SELECT EIGENMODES, GENERATE
1,36,1
** 
*CLOAD, LOAD CASE=1, AMPLITUDE=Powerunit-3376-2_R
_PickedSet32,2,1.
*CLOAD, LOAD CASE=2, AMPLITUDE=Powerunit-3376-2_I
_PickedSet32,2,1.
** ------ Output to data base:
*OUTPUT, FIELD, FREQUENCY=1
*NODE OUTPUT
**V,
UT, 
*END STEP
** 
**
** SIGN: dJrWPqrmhbAQ+Q==
**

随机振动的结果，仿真的位移值在阻尼比0.02时与计算（阻尼比0.16）接近，仿真阻尼比采用0.16时位移的值减小的很大！见下图，原因不明，待查！！！？？？

经过检查，因为仿真计算时激励力的自谱是2，Abaqus仿真用的是1，造成仿真的位移自谱输出是计算值的1/2，另外，计算绘图频率值用的太少，离共振点偏的较大，在阻尼比设成0.02时显示的不对，不是计算有问题，在检查计算公式过程中还发现，结构阻尼系统和粘性阻尼系统的频响函数稍有不同，物体上加激励力的系统与支座上加激励力的频响函数稍有不同，在推导系统响应时，用的相位角phi是滞后角，所以与频响函数中的相位角正负相反，需要注意！！！

随机振动计算的inp文件：

*Heading
** Job name: Job-1 Model name: Model-1
** Generated by: Abaqus/CAE 2016
*Preprint, echo=NO, model=NO, history=NO, contact=NO
**
** PARTS
**
*Part, name=Part-1
*Node
      1,           0.,         100.,           0.
*Nset, nset=Part-1-RefPt_, internal
1, 
*End Part
**  
*Part, name=Part-2
*Node
      1,           0.,           0.,           0.
*Nset, nset=Part-2-RefPt_, internal
1, 
*End Part
**  
**
** ASSEMBLY
**
*Assembly, name=Assembly
**  
*Instance, name=Part-1-1, part=Part-1
*End Instance
**  
*Instance, name=Part-2-1, part=Part-2
*End Instance
**  
*Nset, nset=Set-1, instance=Part-1-1
 1,
*Nset, nset=_PickedSet29, internal, instance=Part-2-1
 1,
*Nset, nset=_PickedSet32, internal, instance=Part-1-1
 1,
*Element, type=MASS, elset=_PickedSet31_Inertia-1_
1, Part-1-1.1
*Mass, elset=_PickedSet31_Inertia-1_
1., 
*Spring, elset=Springs/Dashpots-1-spring

39.4384
*Element, type=SpringA, elset=Springs/Dashpots-1-spring
2, Part-1-1.1, Part-2-1.1
*End Assembly
** ----------------------------------------------------------------
** 
** Random Force Excitation
** 
*Psd-definition, name=Amp-1, TYPE=FORCE
1., 0., 0.1
1., 0., 100.
** 
** BOUNDARY CONDITIONS
** 
** Name: BC-1 Type: Symmetry/Antisymmetry/Encastre
*Boundary
_PickedSet29, ENCASTRE
** Name: BC-2 Type: Symmetry/Antisymmetry/Encastre
*Boundary
**_PickedSet32, YASYMM
_PickedSet32,1,1,0
_PickedSet32,3,6,0
** ----------------------------------------------------------------
** 
** STEP: Step-1-modal
** 
*Step, name=Step-1-modal, nlgeom=NO, perturbation
*Frequency, eigensolver=Lanczos, sim=NO, acoustic coupling=on, normalization=displacement
, , 1000., , , 
** 
** OUTPUT REQUESTS
** 
*Restart, write, frequency=0
** 
** FIELD OUTPUT: F-Output-1
** 
*Output, field
*Node Output, nset=Set-1
UT, 
*End Step
** ----------------------------------------------------------------
** 
** STEP: Step-3-Random Response
** 
*Step, name=Step-3, nlgeom=NO, perturbation
*Random Response
0.1, 2., 20, 3., 
**阻尼由0.02改为0.16,和习题一致
*Modal Damping
1, 2, 0.16
*SELECT EIGENMODES, GENERATE
1,36,1
** 
** BOUNDARY CONDITIONS
** 
** Name: BC-3 Type: Acceleration base motion
**Base Motion, type=ACCELERATION, dof=2, scale=1., load case=3
**Correlation, type=UNCORRELATED, psd=Amp-1, complex=YES
**3, 1., 2.
** 
*Correlation, type=UNCORRELATED, psd=Amp-1, complex=NO
3, 1.
*Cload,load case=3
_PickedSet32,2,1.
** ------ Output to data base:
*OUTPUT, FIELD, FREQUENCY=1
*NODE OUTPUT
**V,
UT, 
*END STEP
** 
**
** SIGN: dJrWPqrmhbAQ+Q==
**