随机分布单向长纤维-复合材料RVE模型-abaqus-python二次开发（上）

前言

之前按照视频教程与书籍写了一篇小学习笔记：Abaqus-python-二次开发学习笔记-单向长纤维RVE，其中只画出了一个纤维，而论文当中，往往是大量随机分布的纤维构成RVE模型，如图所示，其更能反应真实情况。本笔记对此进行了补充学习。

图出自论文
The effects of fiber radius and fiber shape deviations and of
matrix void content on the strengths and failure mechanisms of UD
composites by computational micromechanics

本篇主要参考教程：

https://www.youtube.com/watch?v=UmTwMQDCjG4，在搜寻教程时发现，大部分复合材料RVE教程都是由该博主发布，他应该这方面的专家了，其实他有写好的程序可以直接用，但是是40美元付费软件，目前只能望而却步了，还是需要自己的努力去实现。

ps.鄙人能力水平均有限，免不了做出幼稚的思考和决定，或认知不到位，望诸位和谐讨论，给予批评。

问题分析

关于纤维分布

观察上图的RVE模型，可以看到有很多“切断”的纤维，其上下、左右分别满足周期性。因此在生成随机纤维时，需要考虑中间的普通随机分布的纤维，也需要考虑边界上的纤维以保证周期分布。
观察到图中四个角落上的纤维都只具有一部分，我在此担心后续划分网格时如果纤维此处横截面较小会对网格有影响，因此代码应给予考虑。针对该问题，我在后续的代码实现中直接避免四个角落生成纤维，四边上的纤维也需要满足被切割部分不至于太小的条件。

上述论文采取的是一种叫做nearest neighbor algorithm (NNA)的算法生成随机纤维，体积含量控制在50%。
该方法放到后续再研究实现，目前而言我先采用简单的方法去做实验。

关于abaqusRVE

上述步骤可以在python或matlab中实现，得到纤维的中心坐标，为列表形式，再依次取出在abaqus中创建纤维即可。
此处应该考虑abaqus的代码实现，先GUI界面进行操作，观察代码后再仿照去改写。

GUI试验

与之前的学习文章重复、不重要的内容会略过

创建部件并装配

1. 先画了基体部件：Part-matrix。尺寸：50x50x5
2. 随便画了几个纤维，Part-fiber。尺寸：r=2，depth=5
3. 如图装配

可以看到需要切割纤维突出的部分
.rpy文件代码如下

from abaqus import *
from abaqusConstants import *
session.Viewport(name='Viewport: 1', origin=(0.0, 0.0), width=250.050506591797, 
    height=124.525001525879)
session.viewports['Viewport: 1'].makeCurrent()
session.viewports['Viewport: 1'].maximize()
from caeModules import *
from driverUtils import executeOnCaeStartup
executeOnCaeStartup()
session.viewports['Viewport: 1'].partDisplay.geometryOptions.setValues(
    referenceRepresentation=ON)
# 以上代码为打开abaqus就有的，以下是上述操作生成的
# 首先生成了一个草图'__profile__'定义为s，将草图的几何、顶点、维度、约束分别给到g, v, d, c
s = mdb.models['Model-1'].ConstrainedSketch(name='__profile__', 
    sheetSize=200.0)
g, v, d, c = s.geometry, s.vertices, s.dimensions, s.constraints
# 使ConstrainedSketch对象成为当前视口中的主要对象
s.setPrimaryObject(option=STANDALONE)
# 画一个矩形
s.rectangle(point1=(0.0, 0.0), point2=(50.0, 50.0))
# 新建一个part，是对草图s的拉伸制作，设定深度为5.0
p = mdb.models['Model-1'].Part(name='matrix', dimensionality=THREE_D, 
    type=DEFORMABLE_BODY)
p = mdb.models['Model-1'].parts['matrix']
p.BaseSolidExtrude(sketch=s, depth=5.0)
# 取消ConstrainedSketch对象成为当前视口中的主要对象
s.unsetPrimaryObject()
p = mdb.models['Model-1'].parts['matrix']
session.viewports['Viewport: 1'].setValues(displayedObject=p)
# 草图使用完毕，删除草图
del mdb.models['Model-1'].sketches['__profile__']

# 画纤维，和基体类似
s1 = mdb.models['Model-1'].ConstrainedSketch(name='__profile__', 
    sheetSize=200.0)
g, v, d, c = s1.geometry, s1.vertices, s1.dimensions, s1.constraints
s1.setPrimaryObject(option=STANDALONE)
# 随机纤维的实现将在此处进行更改纤维的坐标
s1.CircleByCenterPerimeter(center=(0.0, 2.0), point1=(2.0, 2.0))
s1.CircleByCenterPerimeter(center=(50.0, 2.0), point1=(52.0, 2.0))
s1.CircleByCenterPerimeter(center=(0.0, 25.0), point1=(2.0, 25.0))
s1.CircleByCenterPerimeter(center=(25.0, 25.0), point1=(27.0, 25.0))
s1.CircleByCenterPerimeter(center=(50.0, 25.0), point1=(52.0, 25.0))
p = mdb.models['Model-1'].Part(name='fiber', dimensionality=THREE_D, 
    type=DEFORMABLE_BODY)
p = mdb.models['Model-1'].parts['fiber']
p.BaseSolidExtrude(sketch=s1, depth=5.0)
s1.unsetPrimaryObject()
p = mdb.models['Model-1'].parts['fiber']
session.viewports['Viewport: 1'].setValues(displayedObject=p)
del mdb.models['Model-1'].sketches['__profile__']

# 装配，定义a为Assembly对象
a = mdb.models['Model-1'].rootAssembly
session.viewports['Viewport: 1'].setValues(displayedObject=a)
session.viewports['Viewport: 1'].assemblyDisplay.setValues(
    optimizationTasks=OFF, geometricRestrictions=OFF, stopConditions=OFF)
a = mdb.models['Model-1'].rootAssembly
a.DatumCsysByDefault(CARTESIAN)
# 将前两个part作为组装体的实例对象
p = mdb.models['Model-1'].parts['fiber']
a.Instance(name='fiber-1', part=p, dependent=OFF)
p = mdb.models['Model-1'].parts['matrix']
a.Instance(name='matrix-1', part=p, dependent=OFF)

根据以上脚本的学习，可以做出简化。

1. 关于session的为视图操作，不必要，全部移除
2. p = mdb.models['Model-1'].parts['fiber']这种多次出现，给予删除

可以看到简化后的代码体量小了很多，在之后的代码中，将直接展示删减后的内容

from abaqus import *
from abaqusConstants import *
from caeModules import *
from driverUtils import executeOnCaeStartup
executeOnCaeStartup()

s = mdb.models['Model-1'].ConstrainedSketch(name='__profile__', 
    sheetSize=200.0)
s.rectangle(point1=(0.0, 0.0), point2=(50.0, 50.0))
p_m = mdb.models['Model-1'].Part(name='matrix', dimensionality=THREE_D, 
    type=DEFORMABLE_BODY)
p_m.BaseSolidExtrude(sketch=s, depth=5.0)
del mdb.models['Model-1'].sketches['__profile__']

s1 = mdb.models['Model-1'].ConstrainedSketch(name='__profile__', 
    sheetSize=200.0)
s1.CircleByCenterPerimeter(center=(0.0, 2.0), point1=(2.0, 2.0))
s1.CircleByCenterPerimeter(center=(50.0, 2.0), point1=(52.0, 2.0))
s1.CircleByCenterPerimeter(center=(0.0, 25.0), point1=(2.0, 25.0))
s1.CircleByCenterPerimeter(center=(25.0, 25.0), point1=(27.0, 25.0))
s1.CircleByCenterPerimeter(center=(50.0, 25.0), point1=(52.0, 25.0))
p_f = mdb.models['Model-1'].Part(name='fiber', dimensionality=THREE_D, 
    type=DEFORMABLE_BODY)
p_f.BaseSolidExtrude(sketch=s1, depth=5.0)
del mdb.models['Model-1'].sketches['__profile__']

a = mdb.models['Model-1'].rootAssembly
a.DatumCsysByDefault(CARTESIAN)
a.Instance(name='fiber-1', part=p_m, dependent=OFF)
a.Instance(name='matrix-1', part=p_f, dependent=OFF)

修剪模型

首先使用合并/切割实体，将之前的两个实例（instance）合并成一个实例，命名为UDcomp。

a.InstanceFromBooleanMerge(name='UDcomp', instances=(a.instances['fiber-1'], 
    a.instances['matrix-1'], ), keepIntersections=ON, 
    originalInstances=SUPPRESS, domain=GEOMETRY)

模型下会自动生成一个部件，接下来双击UDcomp部件，使用切削命令，画两个矩形，进行拉伸切削。


切削后的结果：

# 选取合并体UDcomp的前面和右边，做切割草图的平面。
p = mdb.models['Model-1'].parts['UDcomp']
f, e = p.faces, p.edges
t = p.MakeSketchTransform(sketchPlane=f[23], sketchUpEdge=e[44], 
    sketchPlaneSide=SIDE1, sketchOrientation=RIGHT, origin=(25.0, 25, 
    5.0))
# 草图s按照上述转换生成位置
s = mdb.models['Model-1'].ConstrainedSketch(name='__profile__', 
    sheetSize=141.42, gridSpacing=3.53, transform=t)
# 下面这句根据猜测，为部件根据s草图投影
p.projectReferencesOntoSketch(sketch=s, filter=COPLANAR_EDGES)
# 画一大一小两个矩形，作为切割使用
s.rectangle(point1=(-25.0, -25.0), point2=(25.0, 25.0))
s.rectangle(point1=(-35, -35), point2=(35, 35))
p.CutExtrude(sketchPlane=f[23], sketchUpEdge=e[44], sketchPlaneSide=SIDE1, 
    sketchOrientation=RIGHT, sketch=s, flipExtrudeDirection=OFF)
del mdb.models['Model-1'].sketches['__profile__']

创建材料及设定截面

比较基础，GUI步骤不做展示

mdb.models['Model-1'].Material(name='matrix')
mdb.models['Model-1'].materials['matrix'].Elastic(table=((50.0, 0.3), ))
mdb.models['Model-1'].Material(name='fiber')
mdb.models['Model-1'].materials['fiber'].Elastic(table=((273.0, 0.2), ))
mdb.models['Model-1'].HomogeneousSolidSection(name='Section-matrix', 
    material='matrix', thickness=None)
mdb.models['Model-1'].HomogeneousSolidSection(name='Section-fiber', 
    material='fiber', thickness=None)
p = mdb.models['Model-1'].parts['UDcomp']
c = p.cells
cells = c.getSequenceFromMask(mask=('[#20 ]', ), )
region = regionToolset.Region(cells=cells)
p.SectionAssignment(region=region, sectionName='Section-matrix', offset=0.0, 
    offsetType=MIDDLE_SURFACE, offsetField='', 
    thicknessAssignment=FROM_SECTION)

cells = c.getSequenceFromMask(mask=('[#1f ]', ), )
region = regionToolset.Region(cells=cells)
p.SectionAssignment(region=region, sectionName='Section-fiber', offset=0.0, 
    offsetType=MIDDLE_SURFACE, offsetField='', 
    thicknessAssignment=FROM_SECTION)

mesh

p = mdb.models['Model-1'].parts['UDcomp']
c = p.cells
cells = c.getSequenceFromMask(mask=('[#1f ]', ), )
region = regionToolset.Region(cells=cells)
p.SectionAssignment(region=region, sectionName='Section-fiber', offset=0.0, 
    offsetType=MIDDLE_SURFACE, offsetField='', 
    thicknessAssignment=FROM_SECTION)
p.seedPart(size=1.0, deviationFactor=0.1, minSizeFactor=0.1)

p = mdb.models['Model-1'].parts['UDcomp']
c = p.cells
pickedRegions = c.getSequenceFromMask(mask=('[#3f ]', ), )
p.setMeshControls(regions=pickedRegions, elemShape=WEDGE)
p = mdb.models['Model-1'].parts['UDcomp']
p.generateMesh()

EasyPBC

首先删除两个红叉的实例，否则会报错


最后成功运行

结语

本篇通过简单的案例再次学习了abaqus代码，相比之前新增了切割命令，从而满足周期性条件。
目前存在的问题：指定截面时GUI操作为手动选取，只通过代码无法复现，因此需要新的指定截面的方式，将在下篇更新。