Strain Recorder Problem - Double Tee Girder - Verification purpose

[burakdur]

Hello,

Here is the model code for simply-supported Double Tee girder. Before I move forward to full bridge model, I'd like to verify some fundamental things such as eigenvalues, strain-stress and support reactions etc. This is an elastic analysis, that's why I chose to go with elastic uniaxial material. However, in order to monitor the strain-stress response, I use forceBeamColumn. I haven't defined any pre-stressed tendons or rebars, I'd like to see the outputs for just concrete material first.

I create a fiber section, then, I have top flange + two stems = 3 different patch('quad',.......).

Problems;
1. I can obtain support reaction using recorder, but, why I cannot output them using nodeReaction? I can't also see eleForce outputs.
2. Either using recorder or eleResponse, no way to output stress-strain at the bottom of each stem, the .txt files are all empty and stressStrain gives nothing using eleResponse. Can somebody please explain why this is happening? I tried to make them string and I was playing with it but it did not work out. I also tried to use different beam-column element but they did not work out too.
3. Under a -7.56 N/mm uniform loading on each member (assume that I remove point load at node-6), I am also trying to verify the element forces.

The goal is to obtain the strain-stress response at the last integration of element-5 and first integration of element-6 under a point load of 1000 N at node-6. The fiber location could be right stem bottom mid-point, left stem bottom mid-point or flange top mid-point. You can comment uniform loading case not to run it. The number of fibers are defined separately for the stems and the flange nfZ_top, nfY_top = 5, 20, nfZ_bottom, nfY_bottom = 20, 5.

I'd be glad if anybody help me out to figure out what the problem is. Thank you!

CODE:

def GetModel(nn):
ops.wipe()
ops.model('basic', '-ndm', 2, '-ndf', 3)

#Define unit conversion
inch_to_mm = 25.4
ft_to_mm = 304.8
ksi_to_MPa = 6.89476

Area_Tee = 498 * inch_to_mm * inch_to_mm
L_span = 42 * ft_to_mm

#Define nodes
for i in range(nn+1):
if i == 0:
ops.node(i+1, 0.0, 0.0)
else:
ops.node(i+1, (L_span/nn)*i, 0.0)

#Define support conditions
# tag, DX, DY, RZ
ops.fix(1, 1, 1, 0)
ops.fix(11, 0, 1, 0)

#Calculate uniform distributed loading on the beam
gama = 2400 #kg/m^3
g = 9.80665 #N/kg or m/s^2
Wy = round(gama * g * Area_Tee * (math.pow(10, -9)) * -1, 2) # N/mm

#Calculate the mass
ms = round(gama * Area_Tee * L_span * (math.pow(10, -9)), 2) #kg

#Define nodal masses
for i in range(nn+1):
if i == 0:
ops.mass(i+1, (ms / nn) * 0.5, (ms / nn) * 0.5, (ms / nn) * 0.5)
elif i != 0 and i != nn:
ops.mass(i+1, (ms / nn), (ms / nn), (ms / nn))
else:
ops.mass(i+1, (ms / nn) * 0.5, (ms / nn) * 0.5, (ms / nn) * 0.5)

beamSecTag = 1
transfTag = 1
intTag = 1
matTag_concrete = 1

#Define elastic uniaxial material
#uniaxialMaterial('Elastic', matTag, E, eta=0.0, Eneg=E)
fc = 5000 #psi
Ec = 57 * math.sqrt(5000) * ksi_to_MPa
ops.uniaxialMaterial('Elastic', matTag_concrete, Ec)

nfZ_top, nfY_top = 5, 20
nfZ_bottom, nfY_bottom = 20, 5
y_cg, z_cg = round(30 * inch_to_mm, 2), round(6.87 * inch_to_mm, 2)
tf = round(5 * inch_to_mm, 2)
tw = round(18 * inch_to_mm, 2)

#Flange coordinates in CCW
flange_1_y, flange_1_z = (-1 * y_cg), (-1 * z_cg) + tf
flange_2_y, flange_2_z = (-1 * y_cg), (-1 * z_cg)
flange_3_y, flange_3_z = y_cg, (-1 * z_cg)
flange_4_y, flange_4_z = y_cg, (-1 * z_cg) + tf

dist1 = round(8.75 * inch_to_mm, 2)
dist2 = round(dist1 + 1.0 * inch_to_mm, 2)
dist3 = round(dist2 + 4.50 * inch_to_mm, 2)
dist4 = round(dist3 + 1.0 * inch_to_mm, 2)

#Left web coordinates
Lweb_1_y, Lweb_1_z = y_cg - dist3, tf + tw - z_cg
Lweb_2_y, Lweb_2_z = y_cg - dist4, (-1 * z_cg) + tf
Lweb_3_y, Lweb_3_z = y_cg - dist1, (-1 * z_cg) + tf
Lweb_4_y, Lweb_4_z = y_cg - dist2, tf + tw - z_cg

#Right web coordinates
Rweb_1_y, Rweb_1_z = (-1 * y_cg) + dist2, tf + tw - z_cg
Rweb_2_y, Rweb_2_z = (-1 * y_cg) + dist1, (-1 * z_cg) + tf
Rweb_3_y, Rweb_3_z = (-1 * y_cg) + dist4, (-1 * z_cg) + tf
Rweb_4_y, Rweb_4_z = (-1 * y_cg) + dist3, tf + tw - z_cg

#Define strain-stress fiber coordinates
R_fiby, R_fibz = round(Rweb_1_y + (2.25 * inch_to_mm),2), round(tf + tw - z_cg, 2) #Right stem y and z
L_fiby, L_fibz = round(Lweb_4_y - (2.25 * inch_to_mm),2), round(tf + tw - z_cg, 2) #Left stem y and z
U_flay, U_flaz = 0.0, round(-1 * z_cg,2) #Top flange y and z

ops.section('Fiber', beamSecTag)
#patch('quad', matTag, numSubdivIJ, numSubdivJK, *crdsI, *crdsJ, *crdsK, *crdsL)
ops.patch('quad', matTag_concrete, nfZ_top, nfY_top, flange_1_y, flange_1_z, flange_2_y, flange_2_z, flange_3_y, flange_3_z, flange_4_y, flange_4_z)
ops.patch('quad', matTag_concrete, nfZ_bottom, nfY_bottom, Lweb_1_y, Lweb_1_z, Lweb_2_y, Lweb_2_z, Lweb_3_y, Lweb_3_z, Lweb_4_y, Lweb_4_z)
ops.patch('quad', matTag_concrete, nfZ_bottom, nfY_bottom, Rweb_1_y, Rweb_1_z, Rweb_2_y, Rweb_2_z, Rweb_3_y, Rweb_3_z, Rweb_4_y, Rweb_4_z)

ops.geomTransf('Linear', transfTag)
numIntgrPts = 5

#beamIntegration('Lobatto', tag, secTag, N)
ops.beamIntegration('Lobatto', intTag, beamSecTag, numIntgrPts)

#element('forceBeamColumn', eleTag, *eleNodes, transfTag, integrationTag, '-iter', maxIter=10, tol=1e-12, '-mass', mass=0.0)
eleType = 'forceBeamColumn'
for i in range(nn):
ops.element(eleType, i+1, i+1, i+2, transfTag, intTag, "-mass", 0.0)

ops.recorder('Node', '-file', re[0], '-node', 1,11, '-dof',1,2,3, 'reaction')
ops.recorder('Element', '-file', st[0], '-closeOnWrite','-ele', 5, '-section', str(5), 'fiber', str(R_fiby), str(R_fibz), 'stressStrain')
ops.recorder('Element', '-file', st[1], '-closeOnWrite','-ele', 5, '-section', '5', 'fiber', 'L_fiby', 'L_fibz', 'stressStrain')
ops.recorder('Element', '-file', st[2], '-closeOnWrite','-ele', 6, '-section', '1', 'fiber', 'R_fiby', 'R_fibz', 'stressStrain')
ops.recorder('Element', '-file', st[3], '-closeOnWrite','-ele', 6, '-section', '1', 'fiber', 'L_fiby', 'L_fibz', 'stressStrain')
ops.recorder('Element', '-file', st[4], '-closeOnWrite','-ele', 5, '-section', '5', 'fiber', 'U_flay', 'U_flaz', 'stressStrain')
ops.recorder('Element', '-file', st[5], '-closeOnWrite','-ele', 6, '-section', '1', 'fiber', 'U_flay', 'U_flaz', 'stressStrain')
stressStrain = ops.eleResponse(5, "section", 5, '-fiber', str(-457.2), str(409.7), '1', "stressStrain")

R1 = ops.nodeReaction(1, 2)
R2 = ops.nodeReaction(11, 2)
E1 = ops.eleForce(1)
#Define gravity loading
ops.timeSeries('Linear', 1)
ops.pattern('Plain', 1, 1)
#Define uniform distributed loading
#eleLoad('-ele', *eleTags, '-range', eleTag1, eleTag2, '-type', '-beamUniform', Wy, <Wz>, Wx=0.0, '-beamPoint', Py, <Pz>, xL, Px=0.0, '-beamThermal', *tempPts)
for i in range(nn):
ops.eleLoad('-ele', i+1, '-type', '-beamUniform', Wy, 0.0, 0.0)

#Define nodal load
#ops.load(6, 0.0, -1000.0, 0.0)

ops.constraints('Plain') # how it handles boundary conditions
ops.numberer('Plain') # renumber dof's to minimize band-width (optimization), if you want to
ops.system('BandGeneral') # how to store and solve the system of equations in the analysis
ops.test('NormDispIncr', 1.0e-12, 20)
ops.algorithm('Newton') # use Linear algorithm for linear analysis
ops.integrator('LoadControl', 0.1) # determine the next time step for an analysis, # apply gravity in 10 steps
ops.analysis('Static') # define type of analysis static or transient
ops.analyze(10) # perform gravity analysis

# set damping based on first eigen mode
# -fullGenLapack
eigen = ops.eigen('-genBandArpack', 3)
w, f, T = [], [], []
for k in range(len(eigen)):
w.append(math.pow(eigen[k], 0.5))
f.append(w[k] / (2*math.pi))
T.append(1 / f[k])

values_static = {
"eigen" : eigen,
"Omega" : w,
"frequencies" : f,
"period" : T,
"Mass" : ms,
"Distr. load" : Wy,
"Vertical reaction at left support" : R1,
"Vertical reaction at right support" : R2,
"Element-1 forces": E1,
"Element-5 stress and strain": stressStrain,
"Element y coordinates": R_fiby,
"Element z coordinates": R_fibz
}

for key, value in values_static.items():
print(f"{key} : {value}")

ops.wipe()
#tut +=1

return

[mhscott]

Upgrade to the latest version of OpenSeesPy, then remove all the attempted string conversions from your recorder commands.

[burakdur]

Unfortunately, it didn't solve the issue, still there is nothing being outputted.

[mhscott]

Do the example shown in this post: https://portwooddigital.com/2021/07/25/ ... -response/

[burakdur]

I did that example and it worked, then, I went though my code again and found my mistake while defining the coordinates op.('quad',...), thank you for the help.