This notebook reads files describing a structure, and the files output by Frame2D after an analysis, and checks that the member end forces computed here from the displacements and member loads agree with those computed by Frame2D.
It does this in the simplest way possible, using quite different logic than Frame2D, resulting in a higher degree of confidence in the results. It would have been better had someone else programmed it, but oh well ...
In [31]:
dir = 'KG82'
dir = 'l22x6'
#dir = 'l22x6pd'
def filename(basename):
return dir + '.d/' + basename + '.csv'
def Warn(msg):
print('!!!!! Warning: {}'.format(msg))
In [32]:
import pandas as pd
import math
In [33]:
class Node(object):
def __init__(self,id,x,y):
self.id = id
self.x = x
self.y = y
self.deltaX = 0.
self.deltaY = 0.
self.thetaZ = 0.
In [34]:
ntable = pd.read_csv(filename('nodes'))
NODES = {}
for i,n in ntable.iterrows():
if n.ID in NODES:
Warn("Node '{}' is multiply defined.".format(n.ID))
NODES[n.ID] = Node(n.ID,float(n.X),float(n.Y))
In [35]:
dtable = pd.read_csv(filename('displacements'))
for i,n in dtable.iterrows():
node = NODES[n.ID]
node.deltaX = float(n.DX)
node.deltaY = float(n.DY)
node.thetaZ = float(n.RZ)
In [36]:
pd.DataFrame([vars(v) for v in NODES.values()]).set_index('id')
Out[36]:
In [37]:
dtable
Out[37]:
In [38]:
class Member(object):
E = 200000.
def __init__(self,id,nodej,nodek):
self.id = id
self.nodej = nodej
self.nodek = nodek
dx = nodek.x - nodej.x
dy = nodek.y - nodej.y
self.L = L = math.sqrt(dx*dx + dy*dy)
self.cosx = dx/L
self.cosy = dy/L
self.Ix = 0.
self.A = 0.
self.loads = []
self.releases = set()
for a in 'FXJ FXK FYJ FYK MZJ MZK'.split():
setattr(self,a,0.)
In [39]:
table = pd.read_csv(filename('members'))
MEMBERS = {}
for i,m in table.iterrows():
if m.ID in MEMBERS:
Warn("Member '{}' is multiply defined.".format(m.ID))
MEMBERS[m.ID] = Member(m.ID,NODES[m.NODEJ],NODES[m.NODEK])
In [40]:
import sst
SST = sst.SST()
table = pd.read_csv(filename('properties'))
defIx = defA = None
for i,row in table.iterrows():
if not pd.isnull(row.SIZE):
defIx,defA = SST.section(row.SIZE,'Ix,A')
memb = MEMBERS[row.ID]
memb.Ix = float(defIx if pd.isnull(row.Ix) else row.Ix)
memb.A = float(defA if pd.isnull(row.A) else row.A)
if not pd.isnull(row.Ix):
defIx = row.Ix
if not pd.isnull(row.A):
defA = row.A
In [41]:
table = pd.read_csv(filename('member_loads'))
for i,row in table.iterrows():
memb = MEMBERS[row.ID]
typ = row.TYPE
w1 = None if pd.isnull(row.W1) else float(row.W1)
w2 = None if pd.isnull(row.W2) else float(row.W2)
a = None if pd.isnull(row.A) else float(row.A)
b = None if pd.isnull(row.B) else float(row.B)
c = None if pd.isnull(row.C) else float(row.C)
memb.loads.append((typ,w1,w2,a,b,c))
In [42]:
table = pd.read_csv(filename('releases'))
for i,row in table.iterrows():
memb = MEMBERS[row.ID]
memb.releases.add(row.R.upper())
In [43]:
t = pd.DataFrame([vars(v) for v in MEMBERS.values()]).set_index('id')
del t['nodej']
del t['nodek']
del t['loads']
t
Out[43]:
In [44]:
MEFS = pd.read_csv(filename('mefs')).set_index('ID')
MEFS
Out[44]:
In [45]:
cols = 'FXJ FXK FYJ FYK MZJ MZK'.split()
for m in MEMBERS.values():
for a in cols:
setattr(m,a,0.)
# difference in end displacements, global coords
dX = m.nodek.deltaX - m.nodej.deltaX
dY = m.nodek.deltaY - m.nodej.deltaY
# axial deformation / force:
ldX = dX*m.cosx + dY*m.cosy
T = m.E*m.A*ldX/m.L
m.FXK += T
m.FXJ += -T
#print(m.id,ldX,T)
# shear deformation / force:
vdY = dY*m.cosx - dX*m.cosy
M = -6.*m.E*m.Ix*vdY/(m.L*m.L)
V = 2.*M/m.L
m.MZJ += M
m.MZK += M
m.FYJ += V
m.FYK += -V
#print(m.id,vdY,M,V)
# end rotations / moments:
MJ = (m.E*m.Ix/m.L)*(4.*m.nodej.thetaZ + 2.*m.nodek.thetaZ)
MK = (m.E*m.Ix/m.L)*(2.*m.nodej.thetaZ + 4.*m.nodek.thetaZ)
VJ = (MJ+MK)/m.L
m.MZJ += MJ
m.MZK += MK
m.FYJ += VJ
m.FYK += -VJ
#print(m.id,m.nodej.thetaZ,m.nodek.thetaZ,MJ,MK,VJ)
# applied loads: fixed-end moments and shears:
for ltype,w1,w2,a,b,c in m.loads:
mj = mk = 0.
vj = vk = 0.
if ltype == 'PL':
b = m.L - a
P = w1
mj = -P*a*b*b/(m.L*m.L)
mk = P*b*a*a/(m.L*m.L)
vj = (-P*b + nj + mk)/m.L
vk = -P - vj
elif ltype == 'UDL':
mj = -w1*m.L**2/12.
mk = -mj
vj = -w1*m.L/2.
vk = vj
else:
Warn("Load type '{}' not implemented here ...".format(ltype))
continue
m.MZJ += mj
m.MZK += mk
m.FYJ += vj
m.FYK += vk
# member end moment releases:
relc = m.releases.copy()
if 'MZJ' in m.releases:
mj = -m.MZJ
mk = 0. if 'MZK' in m.releases else 0.5*mj
vj = (mj + mk)/m.L
##print(m.id,'MZJ',m.MZJ,m.MZK,mj,mk,rel)
m.MZJ += mj
m.MZK += mk
m.FYJ += vj
m.FYK += -vj
relc.remove('MZJ')
if 'MZK' in m.releases:
mk = -m.MZK
mj = 0. if 'MZJ' in m.releases else 0.5*mk
vj = (mj + mk)/m.L
##print(m.id,'MZK',m.MZJ,m.MZK,mj,mk,rel)
m.MZJ += mj
m.MZK += mk
m.FYJ += vj
m.FYK += -vj
relc.remove('MZK')
if relc:
Warn("Member end-releases not processed: {}".format(relc))
In [46]:
computed = pd.DataFrame([{k:getattr(m,k) for k in ['id']+cols}
for m in MEMBERS.values()]).set_index('id')
diff = (computed - MEFS[cols])
lim = 1E-12
for c in cols:
biggest = MEFS[c].abs().max()
diff[c][diff[c].abs() < biggest*lim] = 0
diff
Out[46]:
In [47]:
diff.abs().max()
Out[47]:
In [48]:
diff = (computed - MEFS[cols])
for c in cols:
biggest = MEFS[c].abs().max()
r = (diff[c]/biggest)
idr = r.abs().idxmax()
print(c,idr,r[idr])