This program performs an elastic analysis of 2-dimensional structural frames. It has the following features:
In [1]:
from __future__ import print_function
In [2]:
from salib import extend, import_notebooks
import_notebooks()
from Tables import Table
from Nodes import Node
from Members import Member
from LoadSets import LoadSet, LoadCombination
from NodeLoads import makeNodeLoad
from MemberLoads import makeMemberLoad
from collections import OrderedDict, defaultdict
import numpy as np
In [3]:
class Object(object):
pass
class Frame2D(object):
def __init__(self,dsname=None):
self.dsname = dsname
if dsname is not None:
Table.set_source(dsname)
self.reset()
def reset(self):
self.rawdata = Object()
self.nodes = OrderedDict()
self.members = OrderedDict()
self.nodeloads = LoadSet()
self.memberloads = LoadSet()
self.loadcombinations = LoadCombination()
#self.dofdesc = []
#self.nodeloads = defaultdict(list)
#self.membloads = defaultdict(list)
self.ndof = 0
self.nfree = 0
self.ncons = 0
self.R = None
self.D = None
self.PDF = None # P-Delta forces
COLUMNS_xxx = [] # list of column names for table 'xxx'
def get_table(self,tablename,extrasok=False,optional=False):
columns = getattr(self,'COLUMNS_'+tablename)
t = Table(tablename,ds_name=self.dsname,columns=columns,optional=optional)
t.read(optional=optional)
reqdl= columns
reqd = set(reqdl)
prov = set(t.columns)
if reqd-prov:
raise Exception('Columns missing {} for table "{}". Required columns are: {}'\
.format(list(reqd-prov),tablename,reqdl))
if not extrasok:
if prov-reqd:
raise Exception('Extra columns {} for table "{}". Required columns are: {}'\
.format(list(prov-reqd),tablename,reqdl))
return t
In [4]:
%%Table nodes
NODEID,X,Y,Z
A,0.,0.,5000.
B,0,4000,5000
C,8000,4000,5000
D,8000,0,5000
In [5]:
@extend
class Frame2D:
COLUMNS_nodes = ('NODEID','X','Y')
def install_nodes(self):
node_table = self.get_table('nodes')
for ix,r in node_table.data.iterrows():
if r.NODEID in self.nodes:
raise Exception('Multiply defined node: {}'.format(r.NODEID))
n = Node(r.NODEID,r.X,r.Y)
self.nodes[n.id] = n
self.rawdata.nodes = node_table
def get_node(self,id):
try:
return self.nodes[id]
except KeyError:
raise Exception('Node not defined: {}'.format(id))
In [6]:
##test:
f = Frame2D()
In [7]:
##test:
f.install_nodes()
In [8]:
##test:
f.nodes
Out[8]:
In [9]:
##test:
f.get_node('C')
Out[9]:
In [10]:
%%Table supports
NODEID,C0,C1,C2
A,FX,FY,MZ
D,FX,FY
In [11]:
def isnan(x):
if x is None:
return True
try:
return np.isnan(x)
except TypeError:
return False
In [12]:
@extend
class Frame2D:
COLUMNS_supports = ('NODEID','C0','C1','C2')
def install_supports(self):
table = self.get_table('supports')
for ix,row in table.data.iterrows():
node = self.get_node(row.NODEID)
for c in [row.C0,row.C1,row.C2]:
if not isnan(c):
node.add_constraint(c)
self.rawdata.supports = table
In [13]:
##test:
f.install_supports()
In [14]:
##test:
vars(f.get_node('D'))
Out[14]:
In [15]:
%%Table members
MEMBERID,NODEJ,NODEK
AB,A,B
BC,B,C
CD,C,D
In [16]:
@extend
class Frame2D:
COLUMNS_members = ('MEMBERID','NODEJ','NODEK')
def install_members(self):
table = self.get_table('members')
for ix,m in table.data.iterrows():
if m.MEMBERID in self.members:
raise Exception('Multiply defined member: {}'.format(m.MEMBERID))
memb = Member(m.MEMBERID,self.get_node(m.NODEJ),self.get_node(m.NODEK))
self.members[memb.id] = memb
self.rawdata.members = table
def get_member(self,id):
try:
return self.members[id]
except KeyError:
raise Exception('Member not defined: {}'.format(id))
In [17]:
##test:
f.install_members()
f.members
Out[17]:
In [18]:
##test:
m = f.get_member('BC')
m.id, m.L, m.dcx, m.dcy
Out[18]:
In [19]:
%%Table releases
MEMBERID,RELEASE
AB,MZK
CD,MZJ
In [20]:
@extend
class Frame2D:
COLUMNS_releases = ('MEMBERID','RELEASE')
def install_releases(self):
table = self.get_table('releases',optional=True)
for ix,r in table.data.iterrows():
memb = self.get_member(r.MEMBERID)
memb.add_release(r.RELEASE)
self.rawdata.releases = table
In [21]:
##test:
f.install_releases()
In [22]:
##test:
vars(f.get_member('AB'))
Out[22]:
If the SST module is loadable, member properties may be specified by giving steel shape designations (such as 'W310x97') in the member properties data. If the module is not available, you may still give $A$ and $I_x$ directly (it only tries to lookup the properties if these two are not provided).
In [23]:
try:
from sst import SST
__SST = SST()
get_section = __SST.section
except ImportError:
def get_section(dsg,fields):
raise ValueError('Cannot lookup property SIZE because SST is not available. SIZE = {}'.format(dsg))
##return [1.] * len(fields.split(',')) # in case you want to do it that way
In [24]:
%%Table properties
MEMBERID,SIZE,IX,A
BC,W460x106,,
AB,W310x97,,
CD,,
In [25]:
@extend
class Frame2D:
COLUMNS_properties = ('MEMBERID','SIZE','IX','A')
def install_properties(self):
table = self.get_table('properties')
table = self.fill_properties(table)
for ix,row in table.data.iterrows():
memb = self.get_member(row.MEMBERID)
memb.size = row.SIZE
memb.Ix = row.IX
memb.A = row.A
self.rawdata.properties = table
def fill_properties(self,table):
data = table.data
prev = None
for ix,row in data.iterrows():
nf = 0
if type(row.SIZE) in [type(''),type(u'')]:
if isnan(row.IX) or isnan(row.A):
Ix,A = get_section(row.SIZE,'Ix,A')
if isnan(row.IX):
nf += 1
data.loc[ix,'IX'] = Ix
if isnan(row.A):
nf += 1
data.loc[ix,'A'] = A
elif isnan(row.SIZE):
data.loc[ix,'SIZE'] = '' if nf == 0 else prev
prev = data.loc[ix,'SIZE']
table.data = data.fillna(method='ffill')
return table
In [26]:
##test:
f.install_properties()
In [27]:
##test:
vars(f.get_member('CD'))
Out[27]:
In [28]:
%%Table node_loads
LOAD,NODEID,DIRN,F
Wind,B,FX,-200000.
In [29]:
@extend
class Frame2D:
COLUMNS_node_loads = ('LOAD','NODEID','DIRN','F')
def install_node_loads(self):
table = self.get_table('node_loads')
dirns = ['FX','FY','FZ']
for ix,row in table.data.iterrows():
n = self.get_node(row.NODEID)
if row.DIRN not in dirns:
raise ValueError("Invalid node load direction: {} for load {}, node {}; must be one of '{}'"
.format(row.DIRN, row.LOAD, row.NODEID, ', '.join(dirns)))
l = makeNodeLoad({row.DIRN:row.F})
self.nodeloads.append(row.LOAD,n,l)
self.rawdata.node_loads = table
In [30]:
##test:
f.install_node_loads()
In [31]:
##test:
for o,l,fact in f.nodeloads.iterloads('Wind'):
print(o,l,fact,l*fact)
In [32]:
%%Table member_loads
LOAD,MEMBERID,TYPE,W1,W2,A,B,C
Live,BC,UDL,-50,,,,
Live,BC,PL,-200000,,5000
In [33]:
@extend
class Frame2D:
COLUMNS_member_loads = ('LOAD','MEMBERID','TYPE','W1','W2','A','B','C')
def install_member_loads(self):
table = self.get_table('member_loads')
for ix,row in table.data.iterrows():
m = self.get_member(row.MEMBERID)
l = makeMemberLoad(m.L,row)
self.memberloads.append(row.LOAD,m,l)
self.rawdata.member_loads = table
In [34]:
##test:
f.install_member_loads()
In [35]:
##test:
for o,l,fact in f.memberloads.iterloads('Live'):
print(o.id,l,fact,l.fefs()*fact)
In [36]:
%%Table load_combinations
CASE,LOAD,FACTOR
One,Live,1.5
One,Wind,1.75
In [37]:
@extend
class Frame2D:
COLUMNS_load_combinations = ('CASE','LOAD','FACTOR')
def install_load_combinations(self):
table = self.get_table('load_combinations',optional=True)
if len(table) > 0:
for ix,row in table.data.iterrows():
self.loadcombinations.append(row.CASE,row.LOAD,row.FACTOR)
else:
all = self.nodeloads.names.union(self.memberloads.names)
for l in all:
self.loadcombinations.append('all',l,1.0)
self.rawdata.load_combinations = table
In [38]:
##test:
f.install_load_combinations()
In [39]:
##test:
for o,l,fact in f.loadcombinations.iterloads('One',f.nodeloads):
print(o.id,l,fact)
for o,l,fact in f.loadcombinations.iterloads('One',f.memberloads):
print(o.id,l,fact,l.fefs()*fact)
In [40]:
@extend
class Frame2D:
def iter_nodeloads(self,casename):
for o,l,f in self.loadcombinations.iterloads(casename,self.nodeloads):
yield o,l,f
def iter_memberloads(self,casename):
for o,l,f in self.loadcombinations.iterloads(casename,self.memberloads):
yield o,l,f
In [41]:
##test:
for o,l,fact in f.iter_nodeloads('One'):
print(o.id,l,fact)
for o,l,fact in f.iter_memberloads('One'):
print(o.id,l,fact)
In [42]:
##test:
Table.CELLDATA
Out[42]:
In [43]:
@extend
class Frame2D:
def install_all(self):
self.install_nodes()
self.install_supports()
self.install_members()
self.install_releases()
self.install_properties()
self.install_node_loads()
self.install_member_loads()
self.install_load_combinations()
In [44]:
##test:
f.reset()
f.install_all()
In [45]:
@extend
class Frame2D:
def number_dofs(self):
self.ndof = (3*len(self.nodes))
self.ncons = sum([len(node.constraints) for node in self.nodes.values()])
self.nfree = self.ndof - self.ncons
ifree = 0
icons = self.nfree
self.dofdesc = [None] * self.ndof
for node in self.nodes.values():
for dirn,ix in node.DIRECTIONS.items():
if dirn in node.constraints:
n = icons
icons += 1
else:
n = ifree
ifree += 1
node.dofnums[ix] = n
self.dofdesc[n] = (node,dirn)
In [46]:
##test:
f.number_dofs()
In [47]:
##test:
f.ndof, f.ncons, f.nfree
Out[47]:
In [48]:
def prhead(txt,ul='='):
"""Print a heading and underline it."""
print()
print(txt)
if ul:
print(ul*(len(txt)//len(ul)))
print()
In [49]:
@extend
class Frame2D:
def print_nodes(self,precision=0,printdof=False):
prhead('Nodes:')
print('Node X Y Constraints DOF #s')
print('---- ----- ----- ----------- ------')
for nid,node in self.nodes.items():
ct = ','.join(sorted(node.constraints,key=lambda t: Node.DIRECTIONS[t]))
dt = ','.join([str(x) for x in node.dofnums])
print('{:<5s}{:>10.{precision}f}{:>10.{precision}f} {:<11s} {}'\
.format(nid,node.x,node.y,ct,dt,precision=precision))
if not printdof:
return
print()
print('DOF# Node Dirn')
print('---- ---- ----')
for i in range(len(self.dofdesc)):
node,dirn = self.dofdesc[i]
print('{:>4d} {:<4s} {}'.format(i,node.id,dirn))
In [50]:
##test:
f.print_nodes(printdof=True)
In [51]:
@extend
class Frame2D:
def print_members(self,precision=1):
prhead('Members:')
print('Member Node-J Node-K Length dcx dcy Size Ix A Releases')
print('------ ------ ------ ------ ------- ------- -------- -------- ----- --------')
for mid,memb in self.members.items():
nj = memb.nodej
nk = memb.nodek
rt = ','.join(sorted(memb.releases,key=lambda t: Member.RELEASES[t]))
print('{:<7s} {:<6s} {:<6s} {:>8.{precision}f} {:>8.5f} {:>8.5f} {:<10s} {:>10g} {:>10g} {}'\
.format(memb.id,nj.id,nk.id,memb.L,memb.dcx,memb.dcy,str(memb.size),memb.Ix,memb.A,rt,precision=precision))
In [52]:
##test:
f.print_members()
In [53]:
@extend
class Frame2D:
def print_loads(self,precision=0):
prhead('Node Loads:')
if self.nodeloads:
print('Type Node FX FY MZ')
print('---- ---- ---------- ---------- ----------')
for lname,node,load in self.nodeloads:
print('{:<4s} {:<4s} {:>10.{precision}f} {:>10.{precision}f} {:>10.{precision}f}'
.format(lname,node.id,load.fx,load.fy,load.mz,precision=precision))
else:
print(" - - - none - - -")
prhead('Member Loads:')
if self.memberloads:
print('Type Member Load')
print('---- ------ ----------------')
for lname,memb,load in self.memberloads:
print("{:<4s} {:<6s} {}".format(lname,memb.id,load))
else:
print(" - - - none - - -")
prhead("Load Combinations:")
if self.loadcombinations:
print('Case Type Factor')
print('----- ---- ------')
prev = None
for cname,lname,f in self.loadcombinations:
cn = ' '*(len(prev)//2)+'"' if cname == prev else cname
print("{:<5s} {:<4s} {:>6.2f}".format(cn,lname,f))
prev = cname
else:
print(" - - - none - - -")
In [54]:
##test:
f.print_loads()
In [55]:
@extend
class Frame2D:
def print_input(self):
prhead('Frame '+str(self.dsname)+':')
print()
print(' # of nodal degrees of freedom:',self.ndof)
print(' # of constrained nodal degrees of freedom:',self.ncons)
print('# of unconstrained nodal degrees of freedom:',self.nfree,' (= degree of kinematic indeterminacy)')
m = len(self.members)
r = self.ncons
j = len(self.nodes)
c = len(self.rawdata.releases)
print()
print(' # of members:',m)
print(' # of reactions:',r)
print(' # of nodes:',j)
print(' # of conditions:',c)
print(' degree of static indeterminacy:',(3*m+r)-(3*j+c))
print('\n')
self.print_nodes()
print('\n')
self.print_members()
print('\n')
self.print_loads()
In [56]:
##test:
f.print_input()
In [57]:
@extend
class Frame2D:
def buildK(self):
K = np.mat(np.zeros((self.ndof,self.ndof)))
for memb in self.members.values():
Kl = memb.localK()
Tm = memb.transform()
Kg = Tm.T * Kl * Tm
dofnums = memb.nodej.dofnums + memb.nodek.dofnums
K[np.ix_(dofnums,dofnums)] += Kg
return K
In [58]:
##test:
K = f.buildK()
K
Out[58]:
In [59]:
@extend
class Frame2D:
def buildP(self,loadcase='all'):
P = np.mat(np.zeros((self.ndof,1)))
for node,load,factor in self.iter_nodeloads(loadcase):
P[node.dofnums] += load.forces * factor
return P
In [60]:
##test:
f.buildP('one')
Out[60]:
In [61]:
@extend
class Frame2D:
def buildP(self,loadcase='all'):
P = np.mat(np.zeros((self.ndof,1)))
for node,load,factor in self.iter_nodeloads(loadcase):
P[node.dofnums] += load.forces * factor
for memb,load,factor in self.iter_memberloads(loadcase):
fefs = memb.fefs([load])
gfefs = memb.Tm.T * (fefs.fefs * factor)
dofnums = memb.nodej.dofnums + memb.nodek.dofnums
P[dofnums] -= gfefs
return P
In [62]:
##test:
P = f.buildP('one')
P
Out[62]:
In [63]:
@extend
class Frame2D:
def solve(self,loadcase='all'):
self.number_dofs()
K = self.buildK()
P = self.buildP(loadcase)
D = np.mat(np.zeros((self.ndof,1)))
N = self.nfree
Kff = K[:N,:N] # partition the matrices ....
Kfc = K[:N,N:]
Kcf = K[N:,:N]
Kcc = K[N:,N:]
D[:N] = np.linalg.solve(Kff,P[:N] - Kfc*D[N:]) # displacements
R = Kcf*D[:N] + Kcc*D[N:] - P[N:] # reactions at the constrained DOFs
return D,R
In [72]:
##test:
f.reset()
f.install_all()
D,R = f.solve('one')
f.print_input()
In [68]:
##test:
D
Out[68]:
In [69]:
##test:
R
Out[69]:
In [74]:
##test:
# here are the correct reactions:
CR = np.matrix([350E3,412.5E3,-1400E6,0,487.5E3]).T
R-CR
Out[74]:
Comparison is good.
In [ ]: