Statically Determinate 2D Trusses
This notebook developes a function that will compute the forces and reactions in a statically determinate 2D truss and display the results. Inputs are:
- a list of joint names (ids) and x, y coordinates.
- a list of member incidences - the pairs of joints that are connected by truss members.
- a description of the reaction components, giving joint id and direction for each.
- a set of loads expressed as forces on joints, ghiving joint id, magnitude and direction for each.
import numpy as np
from numpy.linalg import solve, matrix_rank
from math import sqrt, atan2, degrees
class SDTError(Exception):
pass
Joints
Each joint object remebers its id and x and y coordinates. A number from 0 to number of joints-1 is sequentially assigned as each joint is created. Each joint is placed in the ‘Joints’ dictionary, indexed by joint id. We check to ensure that all ids are unique, and that no two joints have the same coordinates.
Joints = dict() # global dictionary of Joints
class Joint:
def __init__(self,jid,x,y):
if jid in Joints:
raise SDTError("Joint '{}' is defined more than once.".format(jid))
for oj in Joints.values():
if oj.x == x and oj.y == y:
raise SDTError("Joint '{}' has the same coordinates as '{}'.".format(jid,oj.jid))
self.jid = jid
self.n = len(Joints) # assign sequential number
self.x = x
self.y = y
Joints[self.jid] = self # add new joint to the Joints dictionary
Members
Each truss member is defined by simply naming the joints at either end of the member. An internal id is generated by joing the two joint ids together (i.e., a member from joint d to joint g would have the member id dg. An number from 0 to n-1 is sequentially assigned to each new member in the order they are encountered. Member objects are placed in the global ‘Unknowns’ dictionary, indexed by member id.
Each member has attributes:
- uid - the member id
- n - member (or unknown) number, starting at 0
- ji - the Joint object at the i end (the first end of the pair of joints)
- jj - the Joint object at the j end
- L - the length of the member
- l - a directional cosine - the cosine of the angle between the member axis and the x-axis.
- m - a directional cosine - the cosine of the angle between the member axis and the y-axis.
The following errors are checked as each new Member object is created:
- the joints on either end must be different
- the joint ids have been previously input (Joint objects exist for both ends)
- the member has not been mentioned more than once
- the length of the member is not 0.
Unknowns = dict() # dictionary of objects representing unknowns (member forces and reactions)
JointPairs = dict() # dictionary of pairs of joint ids for asll members - indexed by a tuple of joint ids (jid1,jid2) jid1 <= jid2
class Member:
def __init__(self,jid_i,jid_j):
mid = str(jid_i) + str(jid_j)
if jid_i == jid_j:
raise SDTError("Member '{}' has both joints the same: '{}'.".format(mid,jid_i))
for jid in [jid_i,jid_j]:
if jid not in Joints:
raise SDTError("Joint '{}' for member '{}' is not defined.".format(jid,mid))
if mid in Unknowns:
raise SDTError("Member '{}' is defined more than once.".format(mid))
# remember each pair of connected joints, but make sure the names are in alphabetic order
jpair = (jid_i,jid_j) if jid_i <= jid_j else (jid_j,jid_i)
if jpair in JointPairs:
raise SDTError("Joints '{}' and '{}' are also connected by member '{}'.".format(jid_i,jid_j,JointPairs[jpair].uid))
ji = Joints[jid_i] # retrieve joint objects for both ends
jj = Joints[jid_j]
dx = jj.x - ji.x # coordinate differences
dy = jj.y - ji.y
L = sqrt(dx*dx + dy*dy)
if L == 0.:
raise SDTError("Length of member '{}' is zero.".format(mid))
l = dx/L # direction cosines
m = dy/L
self.uid = mid # store id of member (or 'unknown')
self.n = len(Unknowns) # assign sequential number
self.ji = ji # remeber end joint objects
self.jj = jj
self.L = L # remember geometry
self.l = l
self.m = m
Unknowns[self.uid] = self # store in Unknowns dictionary
JointPairs[jpair] = self # and also in joint pairs dictionary
Reactions
Reactions are also unknowns. They are defined by specifying a direction of the reaction force at a joint. The direction is specified by giving relative coordinate differences along the line of action of the reaction. A direction of 3,4 means the line of action rises 4 units in y for every 3 units in x, and the cosine of the angle between the l.o.a and the x-axis would be 0.6 and the cosine wrt to the y-axis would be 0.8.
The following errors are checked:
- the named Joint has been previously defined in the joint coordinates input
- the direction is valid (i.e, not 0,0)
- the reaction is not defined more than once
Reaction objects have the following attributes:
- jid - the id of the joint
- uid - the id of the reaction (or unknown) - automatically generated from the joint id and direction
- n - the sequential number of the unknown
- ji - the Joint object corresponding to the jid
- jj - None
- l - directional cosine (x-axis) of the line of action
- m - directional cosine (y-axis) of the line of action
Reaction objects are added to the Unknowns dictionary.
class Reaction:
def __init__(self,jid,dx,dy):
if jid not in Joints:
raise SDTError("Joint '{}' is not defined for reaction.".format(jid))
L = sqrt(dx*dx + dy*dy)
if L == 0.:
raise SDTError("dx,dy invalid for joint reaction @ '{}': '{},{}'.".format(jid,dx,dy))
# generate the reaction id from the joint id and direction
rid = 'R' + str(jid)
if dx != 0 and dy == 0:
rid += 'x'
elif dy != 0 and dx == 0:
rid += 'y'
else:
a = degrees(atan2(dy,dx))
if a < 0.:
a += 360.
rid += '{:.0f}'.format(a)
if rid in Unknowns:
raise SDTError("Reaction '{}' already defined at joint '{}'.".format(rid,jid))
self.jid = jid
self.uid = rid
self.n = len(Unknowns) # sequential unknown number
self.ji = Joints[self.jid] # the Joint object
self.jj = None
self.l = dx/L # directional cosines
self.m = dy/L
Unknowns[self.uid] = self
Putting it all together
Function sdtruss() accepts as arguments:
- jc - a list of joint coordinates. Each joint is specified by a 3-tuple of its id (usually a character string) and its x- and y- coordinates.
- mi - a list of member incidences. Each member is specified by a 2-tuple giving the joint ids of the joints at either end of the membr.
- rf - a list of reaction forces. Each reaction is a 3-tuple giving the joint id and the directions of the line of action of the reaction force.
- jl - a list of joint loads. Each load is given as a 4-tuple: the joint id, the magnitude of the load, and the directions of the line of action of the load (directions given as x- and y- relative coordinate differences).
- print_flag - optional, default
True. Whether or not to display a table of resulting meber and reaction forces.
def sdtruss( jc, mi, rf, jl, print_flag=True ):
# empty all the global dictionaries
global Joints, Unknowns, JointPairs
Joints = dict()
Unknowns = dict()
JointPairs = dict()
connected = set() # also the set of joints that have been connected to members
for jid,x,y in jc: # create all the joints from the joint corrdinates list
Joint(jid,x,y)
for jidi,jidj in mi: # create all the members from the member incidences list
Member(jidi,jidj)
connected.add(jidi) # and keep track of connected joints
connected.add(jidj)
for jid in Joints: # now check to ensure that all joints have been connected
if jid not in connected:
raise SDTError("Joint '{}' is not connected to any member.".format(jid))
for jid,dcx,dcy in rf: # create all the reaction unknowns.
Reaction(jid,dcx,dcy)
N = 2*len(Joints) # now check to ensure that number of unknowns matches number of equations
if len(Unknowns) < N:
raise SDTError("Too few unknowns, truss is unstable: {}+{} < 2*{}".format(len(mi),len(rf),len(jc)))
elif len(Unknowns) > N:
raise SDTError("Too many unknowns, truss is statically indeterminate: {}+{} > 2*{}".format(len(mi),len(rf),len(jc)))
# create the C matrix from the direction cosines of the member and reaction forces
C = np.zeros((N,N),np.float64)
for unk in Unknowns.values():
k = unk.n # the column number is the unknown number
ji = unk.ji # Joint objects give the row numbers
jj = unk.jj
C[ji.n*2,k] = unk.l # place the direction cosines
C[ji.n*2+1,k] = unk.m
if jj: # and for members (not reactions), the dcs at the other end
C[jj.n*2,k] = -unk.l
C[jj.n*2+1,k] = -unk.m
# create the load vector from the joint loads input
P = np.zeros((N,1),np.float64)
for jid,p,dcx,dcy in jl: # for each load
if jid not in Joints:
raise SDTError("Joint '{}' not defined in load input.".format(jid))
L = sqrt(dcx*dcx + dcy*dcy)
if L == 0.:
raise SDTError("Improper dcx,dcy for load on joint '{}': {},{}".format(jid,dcx,dcy))
l = dcx/L # direction cosines of the line of action
m = dcy/L
j = Joints[jid] # the joint provides the equation number (row #)
P[j.n*2,0] = -l*p
P[j.n*2+1,0] = -m*p
# now check for solvabilty and solve for forces, if possible
if matrix_rank(C) < N:
raise SDTError("'C' matrix is rank deficient. Truss is unstable.")
Q = solve(C,P)
# optionally display forces and return them
if not print_flag:
return Q
print('unknown joint-i joint-j force')
print('------- ------- ------- -----')
for unk in Unknowns.values(): # for each unknown
uid = unk.uid
iid = unk.ji.jid
jid = unk.jj.jid if unk.jj else '' # reactions don't have a second joint id
t = Q[unk.n,0] # get the force from the unknown number
if jid: # if its a member force, display magnitude and C or T
s = 'T' if t >= 0 else 'C'
print("{:7s} {:7s} {:7s} {:12.4g} {}".format(uid,iid,jid,abs(t),s))
else: # else its a reaction so display one joint + direction
print("{:7s} {:7s} {:7s} {:12.4g} {} @ {:.3g},{:.3g}".format(uid,iid,jid,t,'',unk.l,unk.m))
return Q
Example Usage
if __name__ == '__main__':
# HW2 Q2 2019
jc = [ ('a', 0, 4), # joint coordinates (jid,x,y)
('b', 3, 4),
('c', 9, 4),
('d', 0, 0),
('e', 6, 0),
('f', 9, 0),
]
mi = [ ('a','b'), # member incidences (i,j) joint ids
('b','c'),
('a','d'),
('b','d'),
('a','f'),
('c','e'),
('c','f'),
('d','e'),
('e','f'),
]
rf = [ ('f',1,0), # reaction forces (jid,dcx,dcy) dcx,dcy are direction components
('f',0,1),
('d',0,1),
]
jl = [ ('b',80,0,-1), # joint loads (jid,P,dcx,dcy),
('e',100,0,-1),
]
sdtruss( jc, mi, rf, jl )