abaqus基于python二次开发——钢结构穹顶建模

1. 模型示意
   本工作旨在建立一个上表面近乎球面的钢结构穹顶。如下图所示，该穹顶由环向梁和径向梁组成。环向梁径向梁上下截面都为工字钢。
   
   环向梁截面如下图所示，环向梁截面有一个倾斜角度，为了使其上表面尽可能与球面贴合。
   
   
   径向梁横截面为不经过旋转的工字形
   

2. 代码讲解

 2.1准备工作

在abaqus中运行，需要在开头引入

from abaqus import *
from abaqusConstants import *
from caeModules import *
from driverUtils import executeOnCaeStartup
import numpy as np
executeOnCaeStartup()

定义球体中心点

# center and radius
center = [0,0,0]

定义环向梁径向梁数量与截面尺寸

n=5#环向梁数
m=10#径向梁数
st=0.085#翼缘
sh=0.35/2#腹板

定义球面半径与穹顶底部半径

R=50.235-sh#球面半径
r=21.000-sh#穹顶底部半径

定义径向梁角度和环向梁旋转角度，决定最终建完的穹顶的角度，这里定义为360度，代表建立一个完整的穹顶

ala=2*np.pi
ala2=360

计算径向梁轨迹的角度

an0=np.arcsin(r/R)

2.2 环向梁
一根一根建立环向梁，环向梁根数已于前文定义

for i in range(n):
    angle = an0 * (i+1)/ n
    center = [0, 0]

    s = mdb.models['Model-1'].ConstrainedSketch(name='__sweep__', sheetSize=2000.0)
    g, v, d, c = s.geometry, s.vertices, s.dimensions, s.constraints
    s.setPrimaryObject(option=STANDALONE)

    r0=R*np.sin(angle)
    # 圆弧，圆心，起点、终点、顺时针(CLOCKWISE)，逆时针（COUNTERCLOCKWISE）
    s.ArcByCenterEnds(center=(0, 0), point1=(0, r0),
                      point2=(r0*np.sin(ala), r0*np.cos(ala)), direction=CLOCKWISE)
    # s.CircleByCenterPerimeter(center=(0,0),point1=(0,20))#trajectory
    # s.unsetPrimaryObject()
    s1 = mdb.models['Model-1'].ConstrainedSketch(name='__profile__',
                                                 sheetSize=2000.0, transform=(-0.067804759839806, -0.99769860907143, 0.0,
                                                                             0.0, 0.0, 1.0, -0.99769860907143,
                                                                             0.067804759839806, 0.0,
                                                                             -3.10862446895044e-15, 20.0, 0.0))
    g1, v1, d1, c1 = s1.geometry, s1.vertices, s1.dimensions, s1.constraints
    s1.setPrimaryObject(option=SUPERIMPOSE)
    s1.ConstructionLine(point1=(-100, 0.0), point2=(100, 0.0))
    s1.ConstructionLine(point1=(0.0, -100), point2=(0.0, 100))

    #revolve
    Ibeamx0=np.array([[-1,1],[-1,1],[0,0]])
    Ibeamy0=np.array([[1,1],[-1,-1],[1,-1]])
    Ibeamx=(Ibeamx0-center[0])*np.cos(angle)-(Ibeamy0-center[1])*np.sin(angle)+center[0]
    Ibeamy=(Ibeamy0-center[1])*np.cos(angle)+(Ibeamx0-center[0])*np.sin(angle)+center[1]

    s1.Line(point1=(st*Ibeamx[0][0],st*Ibeamy[0][0]), point2=(st*Ibeamx[0][1],st*Ibeamy[0][1]))
    s1.Line(point1=(st * Ibeamx[1][0], st * Ibeamy[1][0]), point2=(st * Ibeamx[1][1], st * Ibeamy[1][1]))
    s1.CoincidentConstraint(entity1=v1[2], entity2=g1[2], addUndoState=False)
    s1.CoincidentConstraint(entity1=v1[3], entity2=g1[3], addUndoState=False)
    s1.Line(point1=(sh*Ibeamx[2][0],sh*Ibeamy[2][0]), point2=(sh*Ibeamx[2][1],sh*Ibeamy[2][1]))
    s1.PerpendicularConstraint(entity1=g1[4], entity2=g1[6], addUndoState=False)
    s1.CoincidentConstraint(entity1=v1[4], entity2=g1[4], addUndoState=False)
    s1.EqualDistanceConstraint(entity1=v1[0], entity2=v1[1], midpoint=v1[4],
        addUndoState=False)
    s1.CoincidentConstraint(entity1=v1[5], entity2=g1[5], addUndoState=False)
    s1.EqualDistanceConstraint(entity1=v1[2], entity2=v1[3], midpoint=v1[5],
        addUndoState=False)


    p = mdb.models['Model-1'].Part(name='beam'+str(i), dimensionality=THREE_D,
        type=DEFORMABLE_BODY)
    p = mdb.models['Model-1'].parts['beam'+str(i)]
    p.BaseShellSweep(sketch=s1, path=s)
    s1.unsetPrimaryObject()
    p = mdb.models['Model-1'].parts['beam'+str(i)]
    session.viewports['Viewport: 1'].setValues(displayedObject=p)
    name.append('beam'+str(i))
    del mdb.models['Model-1'].sketches['__profile__']
    del mdb.models['Model-1'].sketches['__sweep__']
print(name)

装配环向梁

for j in range