【CivilFEM案例】用Python进行线性和非线性结构分析

这个案例中主要介绍CivilFEM Python脚步在结构建模、加载、分析与结构后处理各个过程中的应用。

案例分析了由杆件和弹簧组成的简单结构，杆件一端固定，另一端与一线性弹簧连接，这端可以在垂直方向上移动，并受到30磅的集中荷载P，方向如下图，加载随时间变化而变化，最后达到30磅。

newDocument("Dim2D","Static","Structural")

ConfigUnits.System = "BInch"

import matplotlib.pyplot as pl

import matplotlib.animation as animation

# Geometry

createPoint("1",Point2d(0,0))

createPoint("2",Point2d(100,1))

createLine("Truss","1","2")

# Materials and sections

createGenericMat("Truss_mat",Double(1e+07),Double(0.3),Double(0),Double(0))

createGenericSection("Truss_sect","Truss_mat",Double(1),Double(1e-10),Double(1e-10),Double(0))

# Structural elements

createSETruss("Truss","Truss_sect","Truss_sect","Truss")

structuralElem0 = StructuralElementsContainer.Find("Truss")

structuralElem0.MeshTool1D.ParameterMesh.NumberOfSegments = 1

createPointSpring("Spring","ToGround",Double(6),"Truss",Point2d(100,1),False,True)

makeZoom("Fit")

# Mesh

mesh()

# Loads and BCs

createLoadGroup("Load_1")

loadGroups4 = LoadGroupsContainer.Find("Load_1")

addPointLoad([loadGroups4],"Load","Truss",Point2d(100,1),Vector2d(0,-1),Double(5))

createLoadGroup("Load_4")

loadGroups4 = LoadGroupsContainer.Find("Load_4")

addPointLoad([loadGroups4],"Load","Truss",Point2d(100,1),Vector2d(0,-1),Double(8))

createLoadGroup("Load_4_8")

loadGroups4 = LoadGroupsContainer.Find("Load_4_8")

addPointLoad([loadGroups4],"Load","Truss",Point2d(100,1),Vector2d(0,-1),Double(12))

createLoadGroup("Load_6_6")

loadGroups4 = LoadGroupsContainer.Find("Load_6_6")

addPointLoad([loadGroups4],"Load","Truss",Point2d(100,1),Vector2d(0,-1),Double(30))

createBCGroup("BC")

bc5 = BCGroupsContainer.Find("BC")

addPointBC([bc5],"1","Truss",Point2d(0,0),True,True)

addPointBC([bc5],"2","Truss",Point2d(100,1),True,False)

# LCs

createLoadCase("LC_1_0s")

loadCase1 = LoadCasesContainer.Find("LC_1_0s")

addLoadGroupToLoadCase([loadCase1],"Load_1",Double(1))

addBCGroupToLoadCase([loadCase1],"BC")

createLoadCase("LC_4_0s")

loadCase1 = LoadCasesContainer.Find("LC_4_0s")

addLoadGroupToLoadCase([loadCase1],"Load_4",Double(1))

addBCGroupToLoadCase([loadCase1],"BC")

loadCase1.LoadCaseCalculationTime = Double(4.0)

createLoadCase("LC_4_8s")

loadCase1 = LoadCasesContainer.Find("LC_4_8s")

addLoadGroupToLoadCase([loadCase1],"Load_4_8",Double(1))

addBCGroupToLoadCase([loadCase1],"BC")

loadCase1.LoadCaseCalculationTime = Double(4.8)

createLoadCase("LC_6_6s")

loadCase1 = LoadCasesContainer.Find("LC_6_6s")

addLoadGroupToLoadCase([loadCase1],"Load_6_6",Double(1))

addBCGroupToLoadCase([loadCase1],"BC")

loadCase1.LoadCaseCalculationTime = Double(6.6)

solve()

# Obtain nodal results for linear case and store them in a list

D_lin=[0.0]*5

Node=Double(0.0)

linear=getResults(r"LC_1_0s.rcf")

linear.getNodeResult(Node, "UTy", 2)

D_lin[1]=-(Node.getValue())

linear=getResults(r"LC_4_0s.rcf")

linear.getNodeResult(Node, "UTy", 2)

D_lin[2]=-(Node.getValue())

linear=getResults(r"LC_4_8s.rcf")

linear.getNodeResult(Node, "UTy", 2)

D_lin[3]=-(Node.getValue())

linear=getResults(r"LC_6_6s.rcf")

linear.getNodeResult(Node, "UTy", 2)

D_lin[4]=-(Node.getValue())

# Add individual solution controls to the LCs for the non linear case

loadCase1 = LoadCasesContainer.Find("LC_1_0s")

loadCase1.SolutionControl.IsParticularControl = True

loadCase1.SolutionControl.InitialTimeStep = Formula("0.1")

loadCase1.SolutionControl.MinimumTimeStep = Formula("0.1")

loadCase1.SolutionControl.MaximumTimeStep = Formula("0.1")

loadCase1 = LoadCasesContainer.Find("LC_4_0s")

loadCase1.SolutionControl.IsParticularControl = True

loadCase1.SolutionControl.InitialTimeStep = Formula("0.1/3")

loadCase1.SolutionControl.MinimumTimeStep = Formula("0.1/3")

loadCase1.SolutionControl.MaximumTimeStep = Formula("0.1/3")

loadCase1 = LoadCasesContainer.Find("LC_4_8s")

loadCase1.SolutionControl.IsParticularControl = True

loadCase1.SolutionControl.InitialTimeStep = Formula("0.125")

loadCase1.SolutionControl.MinimumTimeStep = Formula("0.125")

loadCase1.SolutionControl.MaximumTimeStep = Formula("0.125")

loadCase1 = LoadCasesContainer.Find("LC_6_6s")

loadCase1.SolutionControl.IsParticularControl = True

loadCase1.SolutionControl.InitialTimeStep = Formula("1/18")

loadCase1.SolutionControl.MinimumTimeStep = Formula("1/18")

loadCase1.SolutionControl.MaximumTimeStep = Formula("1/18")

ConfigSolver.CreateIntermediateNonLinearResults = True

ConfigSolver.LargeDeflections = True

solve()

# Add non-linear results to vector

D_n_lin=[0.0]*67

for i in range(1,11):

non_linear=getResults(r"LC_1_0s_"+str(i)+".rcf")

non_linear.getNodeResult(Node, "UTy", 2)

D_n_lin[i]=-(Node.getValue())

for i in range(11,41):

non_linear=getResults(r"LC_4_0s_"+str(i)+".rcf")

non_linear.getNodeResult(Node, "UTy", 2)

D_n_lin[i]=-(Node.getValue())

for i in range(41,49):

non_linear=getResults(r"LC_4_8s_"+str(i)+".rcf")

non_linear.getNodeResult(Node, "UTy", 2)

D_n_lin[i]=-(Node.getValue())

for i in range(49,67):

non_linear=getResults(r"LC_6_6s_"+str(i)+".rcf")

non_linear.getNodeResult(Node, "UTy", 2)

D_n_lin[i]=-(Node.getValue())

# Plot results

fig=pl.figure(figsize=(30, 10), dpi=72)

ax1 = pl.subplot(1, 3, 1)

F_lin=[0, 5, 8, 12, 30]

ax1.plot(D_lin, F_lin, marker='o', color="blue", linewidth=1.0, linestyle="-", label='Linear')

F_n_lin=[]

for x in range (0,11):

F_n_lin.append(5.0/10.0*x)

for x in range (1,31):

F_n_lin.append(3.0/30.0*x+5.0)

for x in range (1,9):

F_n_lin.append(4.0/8.0*x+8.0)

for x in range (1,19):

F_n_lin.append(18.0/18.0*x+12.0)

# First subplot

ax1.plot(D_n_lin, F_n_lin, marker='.', color="red", linewidth=1.0, linestyle="-", label='Non-linear')

ax1.legend(loc='upper left')

ax1.set_ylabel('Force(lbf)')

ax1.set_xlabel('Abs displacement (in)')

ax1.set_title('Linear VS Non-linear')

ax1.set_ylim(0.0, 30.3)

ax1.set_xlim(0.0, 2.7)

ax1.grid(True)

# Animation of second subplot.

ax2 = pl.subplot(1, 3, 2)

ax2.plot([0,100], [0,1], marker='.', color="blue", linewidth=1.0, linestyle="-", label='No Load')

line, =ax2.plot([0,100], [0,1-D_n_lin[i]], marker='.', color="red", linewidth=1.0, linestyle="-", label='Non-linear deformation')

def animate(i):

line.set_data([0,100], [0,1-D_n_lin[i]])

return line,

ax2.set_title('Truss displacement evolution, non-linear case')

ax2.set_ylim(-1.8, 1.1)

ax2.set_xlim(0.0, 101)

ax2.set_ylabel('Real position (in)')

ax2.set_xlabel('Length of truss (in)')

ax2.grid(True)

line_ani = animation.FuncAnimation(fig, animate, frames=67, interval=5, blit=True)