I think it is wrong.
I would be grateful if you could show me where my mistake is.
The image of the curve is available at the link below.
The codes are as follows;
import openseespy.opensees as ops
import opsvis as opsv
import matplotlib.pyplot as plt
import math as m
# wipe model
ops.wipe()
# create model
ops.model('basic', '-ndm', 2, '-ndf', 3)
# define geometric parameters of beams and columns
LCol = 3.5 # m
LBeam = 7 # m
# define section dimensions
HCol = 0.45 # m
BCol = 0.45 # m
HBeam = 0.6 # m
BBeam = 0.45 # m
# create nodes
ops.node(1, 0.0, 0.0)
ops.node(2, 0.0, LCol)
ops.node(3, LBeam, 0.0)
ops.node(4, LBeam, LCol)
ops.node(5, 2*LBeam, 0.0)
ops.node(6, 2*LBeam, LCol)
ops.node(7, 0.0, 2*LCol)
ops.node(8, LBeam, 2*LCol)
ops.node(9, 2*LBeam, 2*LCol)
# define boundary condition
ops.fix(1, 1, 1, 1)
ops.fix(3, 1, 1, 1)
ops.fix(5, 1, 1, 1)
# define section tags
ColSecTag = 1 # assign a tag number to the column section
BeamSecTag = 2 # assign a tag number to the beam section
# define column reinforcement properties
coverCol = 4e-2 # m (cover thickness)
numBarsCol = 8 # number of longitudinal - reinforcement bars in each side of column section (synmetric)
DiamBarsCol = 22e-3 # diameter of column bars
BarAreaCol = ((m.pi*pow(DiamBarsCol, 2))/4) # area of longitudinal - reinforcement bars
# define beam reinforcement properties
coverBeam = 4e-2 # m (cover thickness)
numBarsBeam = 6 # number of longitudinal - reinforcement bars in each side of column section (synmetric)
DiamBarsBeam = 22e-3 # diameter of column bars
BarAreaBeam = ((m.pi*pow(DiamBarsCol, 2))/4) # area of longitudinal - reinforcement bars
# define material tags
IDconc_u = 1 # concrete material tag
IDconc_c_b = 2 # concrete material tag for beam
IDconc_c_c = 3 # concrete material tag for column
IDreinf = 4 # steel material tag
# define concrete properties for Concrete material model(unconfined concrete)
#ops.uniaxialMaterial('Concrete01', IDconc_u, -20e3, -0.0014, -18.63e3, -0.0035)
#ops.uniaxialMaterial('Concrete02', IDconc_u, -20e3, 0.0014, -18.63e3, 3.428e-3, 0.023, 0, 0)
ops.uniaxialMaterial( 'Concrete04' , IDconc_u , -20e3 , -0.002 , -0.0035 , 285e5)
# confined model for beam
# define concrete properties for Concrete material model(confined concrete)
#ops.uniaxialMaterial('Concrete01', IDconc_c_b, -20.10e3, -0.00141, -19.47e3, -5.51e-3)
#ops.uniaxialMaterial('Concrete02', IDconc_c_b, -20.10e3, 0.00141, -19.47e3, 5.51e-3, 0.0054, 0, 0)
ops.uniaxialMaterial( 'Concrete04' , IDconc_c_b , -20.10e3 , -0.00205 , -0.0055 , 285e5)
# confined model for column
# define concrete properties for Concrete01 material model(confined concrete)
#ops.uniaxialMaterial('Concrete01', IDconc_c_c, -20.05e3, -0.00141 , -15.9e3, -5.76e-3)
#ops.uniaxialMaterial('Concrete02', IDconc_c_c, -20.05e3, 0.00141, -15.9e3, 5.76e-3, 0.033, 0, 0)
ops.uniaxialMaterial( 'Concrete04' , IDconc_c_c , -20.05e3 , -0.002025 , -0.0057 , 285e5)
# define reinforcing steel properties for Steel01
#Fy = 420e3 # yield strength (kPa)
#E0 = 2e8 # modulus of elasticity (kPa)
#b = 1e-2 # strain hardening ratio
#ops.uniaxialMaterial('Steel01', IDreinf, Fy, E0, b) # create steel material object
ops.uniaxialMaterial('ReinforcingSteel', IDreinf, 420e3, 525e3, 2e8, 1458.33e3, 0.008, 0.08)
################################################################################################################################
################################################################################################################################
#################################---------------DEFINE FIBER SECTION for COLUMN-----------------################################
#################################-------------------SYNMETRIC SECTION---------------------------################################
################################################################################################################################
################################################################################################################################
# define some parameters for columns
coverY = (HCol/2) # The distance from the section z-axis to the edge of the cover concrete -- outer edge of cover concrete
coverZ = (BCol/2) # The distance from the section z-axis to the edge of the cover concrete -- outer edge of cover concrete
coreY = (coverY - coverCol)
coreZ = (coverZ - coverCol)
nfCoreY = 20 # number of fibers in the core patch in the y direction
nfCoreZ = 20 # number of fibers in the core patch in the z direction
# define fiber section
ops.section('Fiber', ColSecTag)
# define the core patch
ops.patch('quad', IDconc_c_c, nfCoreZ, nfCoreY, -coreY, coreZ, -coreY, -coreZ, coreY, -coreZ, coreY, coreZ)
# Define the four cover patches
ops.patch('quad', IDconc_u, 1, 16, -coverY, coverZ, -coverY, coreZ, coverY, coreZ, coverY, coverZ)
ops.patch('quad', IDconc_u, 1, 16, -coverY, -coreZ, -coverY, -coverZ, coverY, -coverZ, coverY, -coreZ)
ops.patch('quad', IDconc_u, 16, 1, -coverY, coreZ, -coverY, -coreZ, -coreY, -coreZ, -coreY, coreZ)
ops.patch('quad', IDconc_u, 16, 1, coreY, coreZ, coreY, -coreZ, coverY, -coreZ, coverY, coreZ)
# Define reinforcing layers
ops.layer('straight', IDreinf, 3, BarAreaCol, coreY, coreZ, coreY, -coreZ)
ops.layer('straight', IDreinf, 3, BarAreaCol, -coreY, coreZ, -coreY, -coreZ)
ops.layer('straight', IDreinf, 2, BarAreaCol, 0, coreZ, 0, -coreZ)
################################################################################################################################
################################################################################################################################
################################---------------PLOT the FIBER SECTION for COLUMN-----------------###############################
#################################-------------------SYNMETRIC SECTION---------------------------################################
################################################################################################################################
################################################################################################################################
fib_sec_1 = [['section', 'Fiber', ColSecTag],
['patch', 'quad', IDconc_c_c, nfCoreZ, nfCoreY, -coreY, coreZ, -coreY, -coreZ, coreY, -coreZ, coreY, coreZ],
['patch', 'quad', IDconc_u, 1, 4, -coverY, coverZ, -coverY, coreZ, coverY, coreZ, coverY, coverZ],
['patch', 'quad', IDconc_u, 1, 4, -coverY, -coreZ, -coverY, -coverZ, coverY, -coverZ, coverY, -coreZ],
['patch', 'quad', IDconc_u, 4, 1, -coverY, coreZ, -coverY, -coreZ, -coreY, -coreZ, -coreY, coreZ],
['patch', 'quad', IDconc_u, 4, 1, coreY, coreZ, coreY, -coreZ, coverY, -coreZ, coverY, coreZ],
['layer', 'straight', IDreinf, 3, BarAreaCol, coreY, coreZ, coreY, -coreZ],
['layer', 'straight', IDreinf, 3, BarAreaCol, -coreY, coreZ, -coreY, -coreZ],
['layer', 'straight', IDreinf, 2, BarAreaCol, 0, coreZ, 0, -coreZ],
]
matcolor = ['r', 'lightgrey', 'gold', 'w', 'w', 'w']
opsv.plot_fiber_section(fib_sec_1, matcolor = matcolor)
plt.axis('equal')
# plt.savefig('fibsec_rc.png')
plt.show()
################################################################################################################################
################################################################################################################################
#################################---------------DEFINE FIBER SECTION for BEAM-------------------################################
#################################-------------------SYNMETRIC SECTION---------------------------################################
################################################################################################################################
################################################################################################################################
# define some parameters for beams
coverY_Beam = (HBeam/2) # The distance from the section z-axis to the edge of the cover concrete -- outer edge of cover concrete
coverZ_Beam = (BBeam/2) # The distance from the section z-axis to the edge of the cover concrete -- outer edge of cover concrete
coreY_Beam = (coverY_Beam - coverBeam)
coreZ_Beam = (coverZ_Beam - coverBeam)
nfCoreY_B = 20 # number of fibers in the core patch in the y direction
nfCoreZ_B = 16 # number of fibers in the core patch in the z direction
# define fiber section
ops.section('Fiber', BeamSecTag)
# define the core patch
ops.patch('quad', IDconc_c_b, nfCoreZ_B, nfCoreY_B, -coreY_Beam, coreZ_Beam, -coreY_Beam, -coreZ_Beam, coreY_Beam, -coreZ_Beam, coreY_Beam, coreZ_Beam)
# Define the four cover patches
ops.patch('quad', IDconc_u, 1, 16, -coreY_Beam, coverZ_Beam, -coreY_Beam, coreZ_Beam, coreY_Beam, coreZ_Beam, coreY_Beam, coverZ_Beam)
ops.patch('quad', IDconc_u, 1, 16, -coreY_Beam, -coreZ_Beam, -coreY_Beam, -coverZ_Beam, coreY_Beam, -coverZ_Beam, coreY_Beam, -coreZ_Beam)
ops.patch('quad', IDconc_u, 12, 1, -coverY_Beam, coverZ_Beam, -coverY_Beam, -coverZ_Beam, -coreY_Beam, -coverZ_Beam, -coreY_Beam, coverZ_Beam)
ops.patch('quad', IDconc_u, 12, 1, coreY_Beam, coverZ_Beam, coreY_Beam, -coverZ_Beam, coverY_Beam, -coverZ_Beam, coverY_Beam, coverZ_Beam)
# Define reinforcing layers
ops.layer('straight', IDreinf, 3, BarAreaBeam, coreY_Beam, coreZ_Beam, coreY_Beam, -coreZ_Beam)
ops.layer('straight', IDreinf, 3, BarAreaBeam, -coreY_Beam, coreZ_Beam, -coreY_Beam, -coreZ_Beam)
################################################################################################################################
################################---------------PLOT the FIBER SECTION for BEAM-----------------###############################
#################################-------------------SYNMETRIC SECTION---------------------------################################
################################################################################################################################
################################################################################################################################
fib_sec_2 = [['section', 'Fiber', BeamSecTag],
['patch', 'quad', IDconc_c_b, nfCoreZ_B, nfCoreY_B, -coreY_Beam, coreZ_Beam, -coreY_Beam, -coreZ_Beam, coreY_Beam, -coreZ_Beam, coreY_Beam, coreZ_Beam],
['patch', 'quad', IDconc_u, 1, 10, -coreY_Beam, coverZ_Beam, -coreY_Beam, coreZ_Beam, coreY_Beam, coreZ_Beam, coreY_Beam, coverZ_Beam],
['patch', 'quad', IDconc_u, 1, 10, -coreY_Beam, -coreZ_Beam, -coreY_Beam, -coverZ_Beam, coreY_Beam, -coverZ_Beam, coreY_Beam, -coreZ_Beam],
['patch', 'quad', IDconc_u, 8, 1, -coverY_Beam, coverZ_Beam, -coverY_Beam, -coverZ_Beam, -coreY_Beam, -coverZ_Beam, -coreY_Beam, coverZ_Beam],
['patch', 'quad', IDconc_u, 8, 1, coreY_Beam, coverZ_Beam, coreY_Beam, -coverZ_Beam, coverY_Beam, -coverZ_Beam, coverY_Beam, coverZ_Beam],
['layer', 'straight', IDreinf, 3, BarAreaBeam, coreY_Beam, coreZ_Beam, coreY_Beam, -coreZ_Beam],
['layer', 'straight', IDreinf, 3, BarAreaBeam, -coreY_Beam, coreZ_Beam, -coreY_Beam, -coreZ_Beam],
]
matcolor = ['r', 'lightgrey', 'gold', 'w', 'w', 'w']
opsv.plot_fiber_section(fib_sec_2, matcolor = matcolor)
plt.axis('equal')
# plt.savefig('fibsec_rc_beam.png')
plt.show()
# define elastic section for columns
ACol = HCol*BCol
ICol_1 = ((BCol*HCol**3)/12)*0.14
ICol_2 = ((BCol*HCol**3)/12)*0.15
ICol_3 = ((BCol*HCol**3)/12)*0.16
HBeam = 0.60 # m
BBeam = 0.45 # m
ABeam = HBeam*BBeam
IBeam = ((BBeam*HBeam**3)/12)*0.13
ops.section('Elastic', 7, 285e5, ACol, ICol_1) # Column type_1
ops.section('Elastic', 8, 285e5, ACol, ICol_2) # Column type_2
ops.section('Elastic', 9, 285e5, ACol, ICol_3) # Column type_3
ops.section('Elastic', 10, 285e5, ABeam, IBeam) # Beam
# define geomTransf
ColTransfTag = 1 # associate a tag to column transformation
BeamTransfTag = 2 # associate a tag to beam transformation
ColTransfType = 'Linear'
BeamTransfType = 'Linear'
ops.geomTransf(ColTransfType, ColTransfTag) # geometric transfer type for columns
ops.geomTransf(BeamTransfType, BeamTransfTag) # geometric transfer type for beams
# create beamIntegration
Lp_Beam = 0.5*HBeam
Lp_Col = 0.5*HCol
ops.beamIntegration('HingeRadau', 1, ColSecTag, Lp_Col, ColSecTag, Lp_Col, 7) # Column type_1
ops.beamIntegration('HingeRadau', 2, ColSecTag, Lp_Col, ColSecTag, Lp_Col,
# Column type_2ops.beamIntegration('HingeRadau', 3, ColSecTag, Lp_Col, ColSecTag, Lp_Col, 9) # Column type_3
ops.beamIntegration('HingeRadau', 4, BeamSecTag, Lp_Beam, BeamSecTag, Lp_Beam, 10) # Beam
# defined element
ops.element('forceBeamColumn', 1, 1, 2, ColTransfTag, 2)
ops.element('forceBeamColumn', 2, 2, 7, ColTransfTag, 1)
ops.element('forceBeamColumn', 3, 3, 4, ColTransfTag, 3)
ops.element('forceBeamColumn', 4, 4, 8, ColTransfTag, 2)
ops.element('forceBeamColumn', 5, 5, 6, ColTransfTag, 2)
ops.element('forceBeamColumn', 6, 6, 9, ColTransfTag, 1)
ops.element('forceBeamColumn', 7, 2, 4, BeamTransfTag, 4)
ops.element('forceBeamColumn', 8, 7, 8, BeamTransfTag, 4)
ops.element('forceBeamColumn', 9, 4, 6, BeamTransfTag, 4)
ops.element('forceBeamColumn', 10, 8, 9, BeamTransfTag, 4)
# define recorders
# displacement of free nodes
ops.recorder('Node', '-file', 'DİSPLACEMENT.out', '-time', '-node', 7, 8, 9, '-dof', 1, 'disp')
ops.recorder('Node', '-file', 'BASE_SHEAR.out', '-time', '-node', 1, 3, 5, '-dof', 1, 'reaction')
# define vertical load
ops.timeSeries('Linear', 1)
ops.pattern('Plain', 1, 1)
ops.eleLoad('-ele', 7, '-type', '-beamUniform', -23, 0, 0, 7)
ops.eleLoad('-ele', 8, '-type', '-beamUniform', -23, 0, 0, 7)
ops.eleLoad('-ele', 9, '-type', '-beamUniform', -23, 0, 0, 7)
ops.eleLoad('-ele', 10, '-type', '-beamUniform', -23, 0, 0, 7)
# define analysis objects
ops.constraints('Plain')
ops.numberer('Plain')
ops.system('BandGen')
ops.test('NormDispIncr', 1e-7, 10)
ops.algorithm('Newton')
Nstep = 10
Dstep = (1/Nstep)
ops.integrator('LoadControl', Dstep)
ops.analysis('Static')
ops.analyze(Nstep)
# Static Pushover Analysis
ops.wipeAnalysis()
ops.loadConst('-time', 0.0)
IDctrlNode = 7
IDctrlDOF = 1
Dmax = 0.5
Dincr = 0.005
ops.timeSeries('Linear', 2)
ops.pattern('Plain', 2, 2)
ops.load(7, 2, 0.0, 0.0)
ops.load(2, 1, 0.0, 0.0)
# define static pushover analysis objects
ops.constraints('Plain')
ops.numberer('Plain')
ops.system('BandGen')
tol = 1e-4
iter = 1900
ops.test('NormDispIncr', tol, iter, 0, 2)
#ops.test('RelativeNormDispIncr', 1e-7, 2000)
#ops.test('RelativeTotalNormDispIncr', tol, iter)
#ops.algorithm('Newton')
ops.algorithm('Newton')
ops.integrator('DisplacementControl', IDctrlNode, IDctrlDOF, Dincr)
ops.analysis('Static')
Nsteps = int(Dmax/Dincr)
ops.analyze(Nsteps)