在局部-全局分析(图6.8)中,对于全局模型中每个Gauss积分点上给定的应变,使用一个RVE计算应力。全局模型用于计算位移和产生的应变,假设材料是均质的。局部模型考虑不均质性,并通过一个RVE来对非均质材料建模。
图6.8 使用RVE的局部-全局分析
式6.13-6.15在6.2节中用于一次施加一个应变分量,目的是得到材料的等效弹性特性。式6.13-6.15对于施加在RVE的应变的一般状态仍然有效,但是必须注意在边或顶点处的周期性边界的说明。式6.13-6.15为9个约束方程,可以施加在RVE的面上所有成对的周期性的点之间,除了边和顶点处的点。
因为每条边属于两个面,在每条边上,可以看到2个位移分量都被用于施加2个约束方程,每个约束方程来自一个面,如同式6.42-6.44给出。然而,对于一个位移分量,只能写一个约束方程。因此,这些边必须分别处理。类似地,因为3个面交汇于一个顶点,来自3个面的3个周期性约束方程的施加则需要每次使用单个位移分量。根据类似于前面给出的推导,可得下面的结果。
四对顶点需要一次对一个分析。每对顶点关于RVE的中心(坐标为(0, 0, 0))对称。得到的约束方程如下
几何(模块:Part)、材料和截面(模块:Property)、装配(模块:Assembly)、分析步(模块:Step)和网格(模块:Mesh)可以使用例6.2中使用的方法创建,但是要创建2D模型,RVE有1个完整的纤维和4个1/4纤维,如图6.4所示。在本例中没有对称边界条件可以使用。相反,周期性边界条件(PBC)式6.13-6.15必须在Load模块中施加。这些周期性边界条件通过使用约束方程实施。本例中没有载荷本身,而是由式6.13-6.15中等式右侧的常数项(由施加的应变产生)对模型产生加载。
为了写出约束方程,需要定义一个参考点,以及主从节点。然后,主从节点的连接通过精确表示式6.13-6.15的约束方程实现。如前面注意的,顶点的方程与面的方程要分别处理。该工作由[5,PBC_2D.py]自动完成。
求解本例的完整脚本在[5,Ex_6.4-CE.py]可得。首先,该脚本创建网格和节点集 合,稍后用于约束方程。脚本可以在Abaqus/CAE中运行,具体操作为:菜单→File,Run script,[Ex_6.4-CE.py],OK。该脚本调用[5,PBC_2D.py]写约束方程,调用[5,srecover2D.py]计算平均应力。结果在表6.4中给出。
表6.4 在RVE中的平均应力和基体上的最大应力(例6.4)
图6.9 应变的云图(在模型坐标系X-Y-Z中的E11,例6.4)
在微尺度获得RVE所用的过程类似地可以用于在中尺度分析层压板。本例中RVE表示层压板。因此,厚度方向应该保持沿厚度自由膨胀。例如,平行于x-y平面的单层板,并且没有施加式6.15,因此,厚度坐标不受约束(见图6.10)。总之,RVE必须包含整个厚度。对于受面内载荷的对称层压板,RVE可以使用对称边界条件用一半厚度定义(见例6.5)。
图6.10 层压的RVE
和
例6.5 计算[0/90/-45/45]S层压板的Gxy,其中,单层板的材料参数为E1 = 139GPa、E2 = 14.5GPa、G12 = G13 = 5.86GPa、G23 = 2.93GPa、v12 = v13 = 0.21、E1 = 139GPa、v23 = 0.38,单层板厚度tk = 1.25mm。
解6.5 作为层压板对称的一个结果,可以使用层压板厚度一半的代表性体积单元RVE,并在z = 0处使用对称边界条件。RVE的模型使用3D实体单元离散。可以使用Abaqus CAE GUI按照下列伪代码建模。
i. 设置工作路径
菜单:File, Set Work Directory,[C:\Simulia\User\Ex_6.5],OK
菜单:File,Save As,[C:\Simulia\User\Ex_6.5\Ex_6.5.cae],OK
ii. 创建部件
模块:Part
菜单:Part,Create
3D,Deformable,Solid,Extrusion,Approx. size [10],Cont
菜单:Add,Line,Rectangle,[0,0],[1,1],X,Done,Depth [5],OK
菜单:Tools,Datum
Type:Plane,Method:Offset from principal plane
XY Plane,Offset [1.25]
XY Plane,Offset [2.5]
XY Plane,Offset [3.75],X,# 关闭Create Datum对话框
菜单:Tools,Partition
Cell,Use datum plane,# 选取:z=1.25的基准面,
Create Partition
# 选取:较大实体,Done,# 选取:z=2.5的基准面,
Create Partition
# 选取:较大实体,Done,# 选取:z=3.75的基准面,
Create Partition
Done,# 关闭Create Partition对话框
iii. 定义材料、截面并赋截面给部件
模块:Property
菜单:Material,Create
Mechanical,Elasticity,Elastic,Type:Engineering Constants
[139000, 14500, 14500, 0.21, 0.21, 0.38, 5860, 5860, 2930],OK
菜单:Assign,Material Orientation
# 选取:Z=4*1.25(自由表面)处的第1层,Done
Use Default Orientation or Other Method
Definition:Coordinate system,
Additional Rotation Direction:Axis 3
Additional Rotation:Angle [0],OK
# 选取:Z=3*1.25处的第2层,Done
Use Default Orientation or Other Method
Definition:Coordinate system,
Additional Rotation Direction:Axis 3
Additional Rotation:Angle [90],OK
# 选取:Z=2*1.25处的第3层,Done
Use Default Orientation or Other Method
Definition:Coordinate system,
Additional Rotation Direction:Axis 3
Additional Rotation:Angle [-45],OK
# 选取:Z=1*1.25处的第4层,Done
Use Default Orientation or Other Method
Definition:Coordinate system,
Additional Rotation Direction:Axis 3
Additional Rotation:Angle [45],OK,Done
菜单:Section,Create
Solid,Homogeneous,Cont,Material:Material-1,OK
菜单:Assign,Section,# 选取所有,Done,OK,Done
iv. 创建装配
模块:Assembly
菜单:Instance,Create,Independent,OK
v. 定义分析步
模块:Step
菜单:Step,Create
Procedure type:Linear perturbation,Static,Cont,OK
菜单:Output,Field Output Requests,Edit,F-Output-1
Output Variables:[S,E,U,IVOL],OK
vi. 添加边界条件和载荷
# 施加简单剪切
模块:Load
菜单:BC,Manager
Create,Name [yfix],Step:Initial,Disp/Rota,Cont
# 选取:y=0处的表面,Done,# 勾选:U1和U2,OK
Create,Name [disp],Step:Step-1,Disp/Rota,Cont
# 选取:Y=1处的表面,Done
# 勾选:U1 [1],# 勾选:U2 [0],OK
Create,Name [zsymm],Step:Initial,Symm/Anti/Enca,Cont
# 选取:Z=0处的表面,Done,ZSYMM,OK
# 关闭BC Manager对话框
vii. 定义约束
应用Tie约束x = 常数表面,以便得到。
模块:Interaction
菜单:Constraint,Create,Tie,Cont
Choose the master type:Surface,# 选取:表面x=0,Done
Choose the slave type:Surface,# 选取:表面x=1,Done
Position Tolerance:Specify distance [1]
# 不勾选:Adjust slave surface initial position # !!!
# 不勾选:Tie rotational DOFs if applicable,OK
viii. 对模型划分网格
模块:Mesh
菜单:Seed,Instance
Approx global size [1.25],OK
# 对于纯剪,每层厚度一个单元足够
菜单:Mesh,Instance,Yes
ix. 求解并可视化结果
模块:Job
菜单:Job,Manager
Create,Cont,OK,Submit,# 当完成后,Results
菜单:File,Save
菜单:File,Run Script,[Ex_6.5-post.py],OK
脚本[5,Ex_6.5-post.py]可以用于计算平均剪应力,以及剪切模量。因为施加的应变等于单位应变,计算得到的平均应力σ12等于C66。因此,在脚本最后一行得到G12 = 21,441MPa。
""" Post-processing Ex. 6.5 """
# 确认工作路径正确
import os
os.chdir(r’C:\Simulia\User\Chapter_6\Ex_6.5’)
# 打开当前作业的输出数据库(ODB)
from visualization import *
odb = openOdb(path=’Job-1.odb’);
myAssembly = odb.rootAssembly;
# 创建临时变量指向frame库,提高读取速度
frameRepository = odb.steps[’Step-1’].frames;
frameS=[];
frameIVOL=[];
# 在层压板方向(全局)创建坐标系
coordSys = odb.rootAssembly.DatumCsysByThreePoints(name=’CSYSLAMINATE’,
coordSysType=CARTESIAN, origin=(0,0,0),
point1=(1.0, 0.0, 0), point2=(0.0, 1.0, 0.0) )
# 将应力从单层板坐标系转换到定义的CSYSLAMINATE层压板坐标系
# stressTrans=odb.steps[’Step-1’].frames[-1].fieldOutputs[’S’]\
# .getTransformedField(datumCsys=coordSys)
stressTrans=frameRepository[-1].fieldOutputs[’S’]\
.getTransformedField(datumCsys=coordSys)
# 将转换后的应力添加到frameS
frameS.insert(0,stressTrans.getSubset(position=INTEGRATION_POINT));
frameIVOL.insert(0,frameRepository[-1].fieldOutputs[’IVOL’]\
.getSubset(position=INTEGRATION_POINT));
Tot_Vol=0; # 总体积
Tot_Stress=0; # 应力和
for II in range(0,len(frameS[-1].values)):
Tot_Vol=Tot_Vol frameIVOL[0].values[II].data;
Tot_Stress=Tot_Stress frameS[0].values[II].data * frameIVOL[0].values[II].data;
# 计算平均值
Avg_Stress=Tot_Stress/Tot_Vol;
print ’Abaqus/Standard Stress Tensor Order:’
# 根据Abaqus Analysis User’s Manual - 1.2.2 Conventions
print ’ 11-22-33-12-13-23’;
print Avg_Stress;
print ’G12=’,Avg_Stress[3]
练习题
问题6.3 使用问题6.1中计算得到的刚度矩阵,计算E1、E2、v12、v23、G12和G23。
问题6.4 使用问题6.2中计算得到的刚度矩阵,计算E1、E2、v12、v23、G12和G23。
问题6.5 展示例6.3中去掉Tie约束后导致的非零法 正应变和G13的错误结果。
问题6.6 用例6.3中的方法计算G12,但是使用对称边界条件,仅离散RVE的1/4。
问题6.7 通过修改例6.3的边界条件和Tie约束,计算G23。
问题6.8 使用例6.4作为指导,计算G23。