abaqus复合材料微观单胞RVE模型（三种材料）的周期性网格划分，并计算E11,E22,E33,G12,G13和G23等效弹性模量和剪切模量。
刚折腾完Abaqus复合材料RVE建模，必须把踩过的坑给大家盘一盘。这次搞的是包含基体、纤维、界面相的三相材料单胞模型，目标是算出六个方向的等效模量。重点说周期性网格和边界条件这两个要命的部分。
先上硬货——模型分区域划网格才是王道。纤维用六面体扫掠网格，基体用四面体自由网格，界面相用扫掠+切分（图1）。注意了，划完网格必须检查相邻区域的节点是否重叠，用这个Python片段快速验证：
mdb.models['Model-1'].rootAssembly surf_nodes = assembly.surfaces['Fiber_Surf'].nodes base_nodes = assembly.surfaces['Matrix_Surf'].nodes for n1, n2 in zip(surf_nodes, base_nodes): if abs(n1.coords[0] - n2.coords[0]) > 1e-6: print('节点不匹配警告！坐标差:', n1.coords - n2.coords)材料属性设置别手软！纤维用engineering constants各向异性参数，基体用各向同性。关键在界面相——建议用cohesive behavior里的traction-separation定律，厚度设成0.001mm防止数值奇异。
abaqus复合材料微观单胞RVE模型（三种材料）的周期性网格划分，并计算E11,E22,E33,G12,G13和G23等效弹性模量和剪切模量。
周期性边界条件才是核心难点。X方向拉伸工况为例，用方程约束把对应面的位移锁死：
# 创建主从面节点对 master_nset = assembly.sets['X_Master'] slave_nset = assembly.sets['X_Slave'] # 施加位移约束方程 for m_node, s_node in zip(master_nset.nodes, slave_nset.nodes): m_dof = (m_node.label-1)*3 + 1 # U1自由度编号 s_dof = (s_node.label-1)*3 + 1 eq_str = f'1.0*U_{m_dof} -1.0*U_{s_dof} = 0' mdb.models['Model-1'].Equation(name=f'EqX_{m_node.label}', equation=eq_str)这里有个骚操作——用节点编号自动生成约束方程，比手动选点高效十倍。注意节点必须预先按相同顺序排列，不然方程会乱套。
算完六个工况后，用这个脚本自动提取等效模量：
# 读取平均应力和应变 from odbAccess import openOdb odb = openOdb('Job-1.odb') frame = odb.steps['Step-1'].frames[-1] S_avg = frame.fieldOutputs['S'].values[0].data # 应力张量 E_avg = frame.fieldOutputs['E'].values[0].data # 应变张量 # 计算E11等效模量 E11 = S_avg[0]/E_avg[0] if abs(E_avg[0])>1e-8 else 0 print(f'E11等效模量：{E11/1e3:.2f} GPa')特别注意应变接近零时要加判断，否则会出现除以零的报错。建议每个模量单独建分析步，避免交叉干扰。
最后给个实测数据参考：碳纤维/环氧树脂体系E11能到180GPa左右，G23通常最低只有4-5GPa。如果结果偏离太大，先查界面相参数是否合理，再检查周期性约束是否漏节点——这两个坑至少耗了我三天时间...
