粘弹性体接触时接触应力的有限元法分析_许政.pdf
第 6 卷 第 4 期2002 年 12月石河子大学学报(自然科学版)Journal of Shihezi University(Natural Science)Vol.6 No.4Dec.2002文章编号:1007-7383(2002)04-0324-03粘弹性体接触时接触应力的有限元法分析许政1,王庆1,付志一2(1 石河子大学水利建筑工程学院,新疆 石河子 832003;2 中国农业大学工程基础科学部,北京 100083)摘要:以粘弹性Maxwell 模型材料的积分松弛型本构关系为基础,进行了本构关系递推,并建立了Maxwell 模型粘弹性体的有限元分析模型,再利用节点对接触模型理论编制了Maxwell模型二维接触问题接触应力计算有限元程序,并验证了其正确性。最后应用该程序,研究了简单 Maxwell 模型轴对称粘弹性体接触时接触面上轴向应力分布及其变化规律。关键词:粘弹性;接触;有限元法;应力中图分类号:O241.82文献标识码:A粘弹性接触问题在国内外研究都较少。因为接触问题是非线性问题,但它既非材料也非几何非线性问题。在接触问题中,边界条件不是在计算开始前就可给出,它们是计算结果。两接触体间接触面的面积与压力分布随外载变化而变化,并与接触体的刚性有关,这是接触问题的特点,也是它的难点 1。另外对于粘弹性材料,其本构关系较复杂,解析解不易得到。所以,对粘弹性接触问题进行数值方法研究有着重要的理论意义和实际意义。在粘弹性材料的热-湿应力计算中,由于位移、应力和应变都是时间的函数,那么,某一时刻 t 的应力计算必须要用到时间段(0,t)上的所有应变历史。但如果保留每一时间步的应变,必然会占去很大的内存空间和时间。所以,在应力计算中有必要进行递推。我们首先进行了粘弹性松弛型增量本构方程的递推,然后利用节点对接触模型理论编制了Maxwell模型二维接触问题接触应力计算有限元程序,并验证了其正确性,最后应用该程序研究了简单Maxwell模型轴对称粘弹性体接触时接触面上轴向应力分布及其变化规律。利用该程序可计算 Maxwell 模型粘弹性体接触时接触区域大小及接触应力分布情况,这在研究粘弹性体接触数值计算问题上有一定的理论意义。1 粘弹性松弛型增量本构方程的递推粘弹性与弹性力学问题的区别仅在于本构关系不同。从已有的资料来看,生物材料的本构关系一般适于用Maxwell 模型来表示。各向同性粘弹性松弛型本构关系需要由 2 个松弛模量表示,例如剪切松弛模量和体积松弛模量。对于生物材料,由于实验条件困难,通常仅测定其拉压松弛模量,同时假定泊松比为常数。对于具有热-湿流变简单性质的广义 Maxwell 粘弹性体,拉压松弛模量的一般形式为 2,3E(t)=NMi=1Eiexp(-tiMT)。(1)式(1)中,NM 为模型阶数,Ei为模型的第 i 阶刚度,MT为材料的热湿迁移因子。对于许多生物材料,常使用山口信吉的模型,式(1)的一般形式为E(M,T,t)=NMi=1iMTexp(-bitMT)。(2)式(2)中,M 是材料的湿度(含水率),是时间的函繱收稿日期:2002-04-27作者简介:许 政(1970-),男,讲师,硕士,从事固体力学研究。数;T 是材料温度(K),是时间的函数;ai和 bi(i=1 NM)为常数。热湿迁移因子对于稻谷颗粒材料,一般形式为 2,3 MT=exp(cM+dT)。(3)式(3)中,下标 M 表示材料中因水分迁移而产生的湿度效应,下标 T 表示材料中因热传导而产生的温度效应,c 和d 为常数。三向应力状态下,对于具有热-湿流变简单性质的广义 Maxwell 粘弹性体,松弛型积分本构关系式为 4t=SDt-E(Mt,Tt,t-)d。(4)式(4)中,SD=S1S2S2S2S1S20S2S2S1S30S3S3。(5)其中对于三维、空间轴对称、平面应变情况有:S1=1-(1+)(1-2),S2=1-(1+)(1-2),S3=12(1+)。对于平面应力情况有:S1=11-2,S2=1-2,S3=12(1+)。SD 的维数也应做相应变化。为了节省内存空间和时间,在应力计算中有必要进行递推。将时间域用 t0=0,t1,tj-1,tj,tcn=t的划分进行离散,用梯形公式计算积分可得应力增量为cn=S E(cn,cn,0)cn+(1+tcntcn-1)E(cn,cn,tcn)-E(cn-1,cn-1,0)cn-1+S cn-2j=1(1+tj+1tj)-E(cn,cn,tcn-tj)-E(cn-1,cn-1,tcn-tj)j。(6)式(6)中 S=SD2。t=tcn时刻应力的总增量由三部分组成:一是弹性应力增量部分;二是只与上一步应变增量有关部分;三是与第 1 步到第 cn-2 步的应变历史有关的部分,计算需要调用从第 1 步到第 cn-2 步的 j和tj以及 Mcn、Tcn、Mcn-1、Tcn-1、tcn-1等大量信息,此部分需应用递推计算。记第3 部分应力增量为 VVcn,通过计算可得 5VVcn=S cn-2j=1(1+tj+1tj)E(cn,cn,tcn-tj)-E(cn-1,cn-1,tcn-tj)j。(7)对(7)式进行简化,第 1 次引入近似 MTcnMTcn-1,并记icn=i1McnTcnexp(-bitcnMTcn)-1Mcn-1Tcn-1,(8)cn-2i=S CN-2J=1hjexp-bi(tcn-1-tj)MTcn-1 j,(9)可得VVcn=NMi=1icncn-2i。(10)式(7)能否应用递推计算,取决于式(10)能否递推计算。对式(10)分离出 j=cn-2 的项,第 2 次引入近似 MTcn-1MTcn-2,可以得到cn-2i=Shcn-2exp(-bitcn-1MTcn1)cn-2+exp(-bi tcn-1MTcn-1)cn-3i(11)此即为 cn-2i的递推公式。代入式(7)便可得 VVcn的计算公式。根据上述粘弹性增量形式本构方程,应用虚功原理,可以推出增量形式的有限元公式 5。依据这些公式并利用节点对接触模型理论,我们编制了粘弹性Maxwell 体接触问题有限元计算程序,并对程序进行了验证。验证算例之一是用于单向应力接触问题,程序计算结果与解析解误差在 10-6以内;验证算例之二是对文献 6 中刚性圆球与Maxwell 模型半平面接触进行计算,程序计算结果与文献 6 中结果吻合 5。325第 4期许政,等:粘弹性体接触时接触应力的有限元法分析2 粘弹性轴对称Maxwell 体接触时接触面上的应力分析 作为上述粘弹性接触应力的有限元程序应用,我们对4 对简单轴对称Maxwell 体进行计算(计算中时间步长取为松弛时间的 1/10),并对结果加以分析,对数据进行无量纲处理,得到同样的结论 5。图1 中 p 为任一时刻t 时接触压力,p0为 t=0 时接触中心的压力,a 为任一时刻t 时接触圆半径,a0为 t=0 时的接触圆半径。图 1 显示随着时间的推移,持续蠕变的效应是改变压力分布,即由最大压力在接触区中央部分的弹性形式变为压力向边缘集中。这一般在超过松弛时间后出现。图 1轴对称体接触时接触面的应力分布3 结语在研究粘弹性体的接触问题时,因为存在着与时间有关的边界条件,故不能使用“对应原理”去进行解析推导 5,而我们的计算结果也证实了这一点。接触力分布随时间变化,并逐渐远离弹性规律。这一点在粘弹性接触问题的研究中是需要充分注意的。参考文献:1 何君毅,林祥都.工程结构非线性问题的数值解法 M.北京:国防工业出版社,1994.2Shinkinchi Yamachchl,Shingo Ymazawa,KaichiroWak-abayashl,et al.Numerical calculation of internal stresses in brownrice kernel during drying process(part 2)J.农业机械学会杂志,1982,44(1):61-68.3Shinkinchi Yamachchi,Shingoymazawa,KaichiroWak-abayashl,et al.Numerical calculation of internal stresses in brownrice kernel during drying process(part 3)J.农业机械学会杂志,1988,50(2):47-54.4 杨挺青.粘弹性力学 M.武汉:华中理工大学出版社,1990.5 许政.粘弹性体接触有限元法及应用研究 R .北京:中国农业大学,2001.6 K L Johnson.接触力学 M.徐秉业.北京:高等教育出版社,1992.Contact Stress Analysis of Viscoelastic Bodies on Contactwith Finite Element MethodXU Zheng1,WANG Qing1,FU Zhi-yi2(1 College of Water Conservancy&Architecture Engineering,Shihezi University,Shihezi 832003,China;2 Department of Basic Sciences,China Agricultural University,Beijing 100083,China)Abstract:Based on the integral relaxation constitutive relation of Maxwell viscoelastic material,constitutive relation wasdeduced.The analysis model of finite element was suggested.The program of finite element was developed to study thestress of two dimensionsMaxwell viscoelastic bodies on contact by using contact model of node to node,moreover the cor-rectness is validated.Finally,the article calculates the Maxwell model axial symmetry viscoelastic and finds the pressuredistribution rules of interface.Key words:viscoelastic;contact;finite element;stress326石河子大学学报(自然科学版)第 6 卷
