力学大会2015集三维拉压不同模量材料非线性分析的

，（重庆大学航空航天学院，重庆三维拉压不同模量材料的数值分析具有较强的非线性，传统迭代求解过程中刚度矩阵的突变甚至奇异容易导致算法出现收敛为改善算法的收敛性，本文基于参变量变分原理，建立了三维拉压不同模量材料的参变量统一本构方程和参数最小势能原理，并结合有限单元法，将原问题转化为标线性互补问题进而用经典的e算法进过程中，该方法无需根据单元的主应力状态来假设并更新材料的弹性矩阵和单元刚度矩阵现出良好的收敛性和稳成功地用来模拟贝壳类生物复合材料微结构间的单边接触行为。经典弹性力学理论中，基于材料均匀性假设所讨论的是具有相同的拉压弹性模量材料的弹性行为。而在实际工程中，如混泥土、岩石等许多材料在拉伸、压缩时表现出不同的力学响应，一般将这类材料称为拉压不同模量材料。与此同时，随着材料科学和制造技术迅速发展，涌现出大量的新型复合材料，且已经被应用到航空、航天、建筑工程等各个领域。近年来，贝壳类生物复合材料的研究也受到广泛关注，它们有着独特的错缝接合的微观结构，兼具高刚度和高断裂韧性，表现出拉压不同模量材料的特性。对于上述这些材料，如果我们仍采用经典弹性理论来对其进析误差，甚至难以正确描述材料实际的力学行为和响应对这些材料进力学分析的可靠性和准确性，必须采用拉压不同模量弹性理论[1。也逐渐涌现，国内学者在这方面进行了大量的研究。，等构造了拉压不同力问题[6]。2004年，等总结了近年来学者的研究现状，该问题的求解难点[7]。1通讯作者：，E-mail致谢：大连理工大学工业装备结构分析国家开放基金No.国家自然科学基金NaturalScienceFoundationProject 问题量材料的小变形算例[9]和大位移算例[10,11]证实了PVP算法确实能够大大提高收敛外，Du和Guo提出了一系列拉压不同模量问题的功能方程，还较为准确的估计出了仿复合纳米材料的弹性边界[13]以获得相应的弹性矩阵，然后进入下一次迭代。其中存在的可归纳为两点：(1)服这两点。因而，发展针对拉压不同模量力学分析的高效、稳定计算方法仍然十分三维拉压不同模量弹性问题的平衡微分方程位移应力分量u是位移分量b为体力L和A为微分算子：Aσb=ε 0 0 0 0

A

0 0

0 0 拉压不同模量材料的应力应变关系通常可以简化为一个双线性模型（如图1所示），材料的拉伸模量为E，压缩模量为E，拉压不同模量材料的应变是不连续的。此外，将三维拉压不同模量材料的本构方程建立在材料的主应力方向上。如图2所示，图中OxIyIzI为主方向的坐标系。主方向的应力和应变和柔度矩阵分别用σI、εI和aI来表示。三维拉压不同模量材料在主应力方向σσεσε 图1双线性本构关系(a)EE(b)EεIaIσI或σI ε ε γσ σ γ

aaD1 a a εi(iα,β,γ)、σi(iα,β,γ)分别表示主方向的应变分量和应力分量。DI为主方向的弹性矩阵。每个单元柔度矩阵aI中的柔度系数由该单元主应力的正负符号来决定。所以，、ν分别为拉压方向的泊松比。此外，基于公共假设νEνE，可以保证柔度矩阵a1,σ

(i a1,σ aijE,σj ν

a ,σ 在推导含参变量的本构方程时，要分别考虑EE和EE两种情况。以下EE进行详细的推导，对于EE的情形，推导过程是类似的。当EE时，可以由料的柔度矩阵a。在三个主应力方向上分别引入参数变量λα、λβ和λγ，所以在八种状态

E

ε

ε ε E

σλγ

β

1γ

E EEEEEE

σ σ σ EEλ E EEλ EE 由于EE，所以EE0。可以确定λαλβ、λγ的值于等于0。事实上，λα的值表示σα/E与σα/E之差。将E假定为σα0图3所示）。因此，λα0xI上处于拉伸状态，λα0则表示材料在该方将松弛变量να，νβ和νγ分别引入方程(11)-(13)后，得到等效的互补方程如下所示：E

σIEEλνλi λ λT，ν ν γ

E DIεIEEIDIλν

0

EI=0 0， DI 00

1EE EE综合EE和EE两种情况，将拉压不同模量材料的统一参变量本构方程写

sDεsδIDλν I λi ssign(EE

EEEE

sign(a)

ifaifa

应变关系。只要参数变量λ，松弛变量ν以及物理变量同时满足方程(19)，方程材料的应变能相对于材料所处的参考系是独立的，所以可以直接在主方向坐标OxIyIzI中表示出材料的应变能，材料的应变 ψ1 2

sλTDεsλ ψ1QεsλTDQεsλ min.λ QLusλTDQLusλdV uTbdV 2 s.t.sDεsδIDλν I λi 可以看出，参数变量λ应力所处的状态对平衡方程的变分计算之中[11(24)(24b)变ν和参数变量λ是互补的。不难证明方程(24)的解即为三维拉压不同模量问题在已知应力和位移边界条件下的真实位移。有关于一般弹塑性问题的参变量变分原理已经得到证明[11压不同模量问题的证明过程与其相似，有的读者可参考相关文献。值求解。离散后的势能表达式为

1qTKqqTWλ2

KKe

WWe

qqe

λλe

f

E 是结构单总数，基于参变量变分原理和有限元平衡方程，对方程(25)进变计可qK1(Wλf HqAλν 0,λNG表示结构中所有点的数目H DQ

I

AdiagsδIDI

λi求解[11]。结合方程(27)和(31)，通过对λ的迭代计算，位移q可以求出，这样一来，可以解决三维拉压不同模量材料的非线性分析问题。当前后两次的迭代误qn

ψ算例1如图2b所示的三维H10.0m，正方形横截面边长W1.0m，顶P=6.0Nm2，材料密度γ2.0N/m3底部受固定约束。用80个八节点实体元来理该该结构，若结构的泊松比为0，则该结构可转换为图2a所示的一维杆.分析以下几种状况：(a)EE5000Pa;(b)E1000Pa,E5000Pa(c)E500Pa,E5000Pa 图2(a)一维杆件;(b)三维柱体;(c)一压分布;(b)三压分画出了杆件的拉压分布图，可见一维和三维时拉压分布是一致的。分界面都在w7.0m

E/E

E/E

E/E 图3门字型结用384个八节点三维实体单元离散该门字型结构。压缩E1800MPa，压缩泊松比ν0.3。分析以下五种情况：(a)E1800MPaE/E1；(b)E900MPa，E/E2；(c)E600MPa，E/E3；(d)E450MPa，E/E4，(e)E360MPaE/E5。表2中列出了点A在拉压模量比不同的五种情况下用三种算法解得的A1-5.462e--5.462e-1-5.462e--5.462e--5.462e-2-6237e---3-6.827e---

表2A4-7329e- 5-7.775e- 当EE时，NR、MNR方法均无法求出结果，PVP方法可以成功算出结果。且随着E/E值的逐渐变大，A点在竖直方向的位移也逐渐增大，又以E/E3为例，绘制出图4各算法的收敛曲线E/E3算更为复杂的结构，抑或是结构的大变形这类有着更加非线性特性的问题时，相算例图5贝壳显微结构图 图6贝壳类生物复合材料的示意图(a)砖墙模型(b)典型代表体积元(c)受拉力作用(d)受压成。无机物的杨氏模量Em100GPa，有机物材料的杨氏模量则被设定为从1到10GPa。材料的拉伸和压缩模量分别为E0.1GPa、E1000GPa，泊松比都为0.33。无机物材料的体积占了所有材料总体积的90%，无机物材料和有机物材料的宽度分别为hm0.5μm、ho20nm，用PVP算法对该材料进行数值模拟，使其进行单轴的拉伸和压缩，之后便可通过ER/来计算出该生物复合材料的等效模量（E表示材料的等效模量，为自由端受位移约束的强迫位移，R为约束端受到的等效外力）。图7包含三种图8拉伸和压缩图9等效模量位于上、下界之在图8中，拉伸和压缩方向的等效模量分别EtEc来表示，下标1 等分析模型求出的等效模量基本一致[14]。此外，该结果也落在由Du和Guo[13]用变AmbartsumyanC.A,著,不同模量弹性理论[M].,译,:中国铁道,TimoshenkoS.,Strengthofmaterials.PartII:advancedtheoryandproblems[M].2nded.Princeton:VanNostrand,1941.3,王志峰,不同拉压模量弹性力学问题的有限元法[J].计算结构力学及其应用,1989,6:236-杨海天 ,杨克俭,等,初应力法解拉压双弹性模量问题[J].大连理工大学学报,1992,32:35-杨海天 ,光滑函数法求解拉压不同模量问题[J].计算力学学报,2006,23:19-YangH.T.,WangB.,Anysisoflongitudinalvibrationofbimodularrodviasmoothingfunctionapproach[J].JournalofSoundandVibration,2008,317:419-431.7,,,不同模量弹性问题理论及有限元法研究进展[J].力学与实践,2004,26:9-HeX.T.,ZhengZ.L.,SunJ.Y.,etal,Convergenceysisofafiniteelementmethodbasedondifferentmoduliintensionandcompression[J].InternationalJournalofSolidsandStructures,2009,46:3734-3740.ZhangL.,GaoQ.,ZhangH.W.,Anefficientalgorithmformechanicalysisofbimodulartrussandtensegritystructures[J].Int.J.Mech.Sci.,2013,70:57-68.ZhangH.W.,ZhangL.,GaoQ.,Numericalmethodfordynamicysisoftwo-dimensionalbimodularstructures[J].AIAAJ.,2012,50:1033-1942.11,张洪武,,参变量变分原理及其在工程中的应用[M].:科学,张洪武,参变量变分原理与材料和结构力学分析[M].:科学,DuZ.L.,GuoX.,Variationalprinciplesandtherelatedboundingtheoremsforbi-modulusmaterials[J].Mech.Phys.Solids,2014,73:183-211.14.BertoldiK.,BigoniD.,DruganW.J.,Nacre:anorthotropicandbimodularelasticmaterial[J].Compos.Sci.Technol,2008,68(6):1363-1375.QuadraticProgrammingAlgorithmforNonlinearysisof3D-BimodularZAHNGHui-ting,ZHANGLiang,YAN(CollegeofAerospaceEngineering,ChongqingUniversity,Chongqing400030, ysisof3Dbimodularmaterialsexhibitsstrongnonlinearityandconvergenceisusuallyaproblemfortraditionali tivesolution.Toimproveconvergenceofthealgorithm,aunifiedconstitutiveequationof3Dbimodularmaterialsandacomplementaryfiniteelementformulationareproposedbasedonparametricvariationalprinciple(PVP),andtheproblemof3D