Cách thực hiện phân tích đẩy bằng Python

lập trình

[ad_1]

Tôi đang thực hiện phân tích đẩy bằng cách sử dụng Openees trong Spyder. Nhưng tôi chỉ nhận được đường cong hiệu suất tuyến tính. Bạn có thể giúp tôi tìm ra vấn đề ở đâu không?

Những gì tôi đã thử:

print("==========================")

from openseespy.opensees import *


import os, sys
import numpy as np
import math as math
import matplotlib.pyplot as plt


print("Starting RCFrameGravity example")

wipe()
# Create ModelBuilder (with two-dimensions and 3 DOF/node)
model('basic', '-ndm', 2, '-ndf', 3)

# Create nodes
# ------------

# Set parameters for overall model geometry
height = 1300.0

# Create nodes
#    tag, X, Y
node(1, 0.0, 0.0)
node(2, 0.0, 0.0)
node(3, 0.0, height)

# Fix supports at base of columns
#   tag, DX, DY, RZ
fix(1, 1, 1, 1)
fix(2, 1, 1, 0)

#%%
# Define materials
uniaxialMaterial('Hysteretic', 1, 0.006, 75514606.25, 0.03, 75514606.25*1.1, 
                 -0.006, -75514606.25, -0.03, -75514606.25*1.1, 1.0, 1.0, 0.0, 0.0, 0.7)
#%%
# Define elements
# ----------------------

# Geometry of column elements
#                tag

geomTransf('PDelta', 1)

# Number of integration points along length of element
np = 10

# Lobatto integratoin
beamIntegration('Lobatto', 1, 1, np)

# Create the coulumns 
element('zeroLength',  1,   1,   2,  '-mat',1, '-dir',6)
bc = 250.0
dc = 250.0
A  = bc*dc
Iy = bc*dc**3/12
Ec = 26692.696
element('elasticBeamColumn',  2,   2,   3,  A,  Ec,  Iy, 1)

# import vfo rendering module
import vfo.vfo as vfo

# render the model after defining all the nodes and elements
# vfo.plot_model(show_eletags="yes")

# plot mode shape
# vfo.plot_modeshape(modenumber=3)

#%%

# Define gravity loads
# --------------------

#  a parameter for the axial load
P = 200000.0

# Create a Plain load pattern with a Linear TimeSeries
timeSeries('Linear', 1)
pattern('Plain', 1, 1)

# Create nodal loads at nodes 3 
#    nd  FX,  FY, MZ
# load(3, -P, 0.0, 0.0)
load(3, 0.0, -P, 0.0)

# ------------------------------
# End of model generation
# ------------------------------


# ------------------------------
# Start of analysis generation
# ------------------------------

# Create the system of equation, a sparse solver with partial pivoting
system('ProfileSPD')

# Create the constraint handler, the transformation method
constraints('Transformation')

# Create the DOF numberer, the reverse Cuthill-McKee algorithm
numberer('RCM')

# Create the convergence test, the norm of the residual with a tolerance of
# 1e-12 and a max number of iterations of 50
test('NormDispIncr', 1.0e-12, 50, 3)

# Create the solution algorithm, a Newton-Raphson algorithm
algorithm('Newton')

# Create the integration scheme, the LoadControl scheme using steps of 0.1
integrator('LoadControl', 0.1)

# Create the analysis object
analysis('Static')

# ------------------------------
# End of analysis generation
# ------------------------------


# ------------------------------
# Finally perform the analysis
# ------------------------------

# perform the gravity load analysis, requires 50 steps to reach the load level
analyze(50)

# Print out the state of nodes 3 
# print node 3 

# Print out the state of element 1 and 2
# print ele 1 2

u3 = nodeDisp(3, 1)

results = open('results.out', 'a+')

results.close()

print("==========================")

#%%

# ----------------------------------------------------
# Start of Model Generation & Initial Gravity Analysis
# ----------------------------------------------------

# Set the gravity loads to be constant & reset the time in the domain
loadConst('-time', 0.0)

# ----------------------------------------------------
# End of Model Generation & Initial Gravity Analysis
# ----------------------------------------------------


# ----------------------------------------------------
# Start of additional modelling for lateral loads
# ----------------------------------------------------

# Define lateral loads
# --------------------

# Set some parameters
H = 200000.0  # Reference lateral load

# Set lateral load pattern with a Linear TimeSeries
pattern('Plain', 2, 1)

# Create nodal loads at nodes 3
#    nd    FX  FY  MZ
load(3, H, 0.0, 0.0)

# ----------------------------------------------------
# End of additional modelling for lateral loads
# ----------------------------------------------------

# ----------------------------------------------------
# Start of modifications to analysis for push over
# ----------------------------------------------------

# Set some parameters
dU = 1.0  # Displacement increment

# Change the integration scheme to be displacement control
#                             node dof init Jd min max
integrator('DisplacementControl', 3, 1, dU, 1, dU, dU)

# ----------------------------------------------------
# End of modifications to analysis for push over
# ----------------------------------------------------

# ------------------------------
# Start of recorder generation
# ------------------------------

# Stop the old recorders by destroying them
# remove recorders

# Create a recorder to monitor nodal displacements
# recorder Node -file node3.out -time -node 3 -dof 1 2 3 disp

# Create a recorder to monitor element forces in columns
# recorder EnvelopeElement -file ele12.out -time -ele 1 2 forces

# --------------------------------
# End of recorder generation
# ---------------------------------

# ------------------------------
# Finally perform the analysis
# ------------------------------

# Set some parameters
maxU = 250.0  # Max displacement
currentDisp = 0.0
ok = 0

test('NormDispIncr', 1.0e-12, 1000)
algorithm('ModifiedNewton', '-initial')


topdispl = [] #empty list to record top node displacement
ndispl = {1: [], 2: [], 3: []}
eforce = {(1,1): [], (1,2):[], (1,3):[], (1,4):[], (1,5):[], (1,6):[], 
          (2,1): [], (2,2):[], (2,3):[], (2,4):[], (2,5):[], (2,6):[]}

while ok == 0 and currentDisp < maxU:

    ok = analyze(1)
    print('currentDisp =',currentDisp)

    # if the analysis fails try initial tangent iteration
    if ok != 0:
        print("modified newton failed")
        break
    # print "regular newton failed .. lets try an initail stiffness for this step"
    # test('NormDispIncr', 1.0e-12,  1000)
    # # algorithm('ModifiedNewton', '-initial')
    # ok = analyze(1)
    # if ok == 0:
    #     print "that worked .. back to regular newton"

    # test('NormDispIncr', 1.0e-12,  10)
    # algorithm('Newton')

    currentDisp = nodeDisp(3, 1)
    
    #record each step
    topdispl.append(nodeDisp(3,1))
    for ntag in [1,2,3]:
        ndispl[ntag].append(nodeDisp(ntag,1))
        
    for etag in [1,2]:
        eforce[etag,1].append(eleForce(etag,1))
        eforce[etag,2].append(eleForce(etag,2))
        eforce[etag,3].append(eleForce(etag,3))
        eforce[etag,4].append(eleForce(etag,4))
        eforce[etag,5].append(eleForce(etag,5))
        eforce[etag,6].append(eleForce(etag,6))
       

xpoints=ndispl[3]
ypoints=eforce[2,1]
plt.plot(xpoints, ypoints)

    
results = open('results.out', 'a+')

results.close()

# Print the state at node 3
# print node 3

print("==========================")

[ad_2]

コメント

タイトルとURLをコピーしました