from AppendixAcode import *
from AppendixDcode import *
class wfsection():
# -------------------------------------------------
# Constructor
# -------------------------------------------------
def __init__(self, hw, bf, tf, tw, nf, nw, material):
# Input
self.mat = material
self.hw = hw
self.bf = bf
self.tf = tf
self.tw = tw
self.nf = nf
self.nw = nw
# Storing fibers and positions
self.fibers = []
self.locations = []
self.areas = []
# Create flange fibers
fibthick = self.tf / self.nf
area = fibthick * self.bf
yf = 0.5*self.hw - 0.5*fibthick
for i in range(0, self.nf*2, 2):
# Top flange
if self.mat[0] == 'Bilinear':
fiber = bilinearMaterial(self.mat)
elif self.mat[0] == 'BoucWen':
fiber = boucWenMaterial(self.mat)
else:
print("Error: Wrong material given to wide-flange section")
import sys
sys.exit()
self.fibers.append(fiber)
self.locations.append(yf)
self.areas.append(area)
# Bottom flange
if self.mat[0] == 'Bilinear':
fiber = bilinearMaterial(self.mat)
elif self.mat[0] == 'BoucWen':
fiber = boucWenMaterial(self.mat)
else:
print("Error: Wrong material given to wide-flange section")
import sys
sys.exit()
self.fibers.append(fiber)
self.locations.append(-yf)
self.areas.append(area)
yf = yf - fibthick
# Create web fibers
fibthick = (self.hw - 2 * self.tf) / self.nw
area = fibthick * self.tw
yf = 0.5*self.hw - self.tf - 0.5*fibthick
for i in range(self.nw):
if self.mat[0] == 'Bilinear':
fiber = bilinearMaterial(self.mat)
elif self.mat[0] == 'BoucWen':
fiber = boucWenMaterial(self.mat)
else:
print("Error: Wrong material given to wide-flange section")
import sys
sys.exit()
self.fibers.append(fiber)
self.locations.append(yf)
self.areas.append(area)
yf = yf - fibthick
# Number of fibers
self.numfib = len(self.fibers)
# Container for section deformations
self.vs = [0.0, 0.0]
# Commit total deformations to accommodate element13
self.trialUsTotal0 = 0.0
self.trialUsTotal1 = 0.0
self.committedUsTotal0 = 0.0
self.committedUsTotal1 = 0.0
# -------------------------------------------------
# Set parameter
# -------------------------------------------------
def setParameter(self, parameter, value):
if parameter == 'hw':
self.hw = value
elif parameter == 'bf':
self.bf = value
elif parameter == 'tf':
self.tf = value
elif parameter == 'tw':
self.tw = value
else:
for i in range(self.numfib):
self.fibers[i].setParameter(parameter, value)
return
# Storing fibers and positions
self.fibers = []
self.locations = []
self.areas = []
# Create flange fibers
fibthick = self.tf / self.nf
area = fibthick * self.bf
yf = 0.5*self.hw - 0.5*fibthick
for i in range(0, self.nf*2, 2):
# Top flange
if self.mat[0] == 'Bilinear':
fiber = bilinearMaterial(self.mat)
elif self.mat[0] == 'BoucWen':
fiber = boucWenMaterial(self.mat)
else:
print("Error: Wrong material given to wide-flange section")
import sys
sys.exit()
self.fibers.append(fiber)
self.locations.append(yf)
self.areas.append(area)
# Bottom flange
if self.mat[0] == 'Bilinear':
fiber = bilinearMaterial(self.mat)
elif self.mat[0] == 'BoucWen':
fiber = boucWenMaterial(self.mat)
else:
print("Error: Wrong material given to wide-flange section")
import sys
sys.exit()
self.fibers.append(fiber)
self.locations.append(-yf)
self.areas.append(area)
yf = yf - fibthick
# Create web fibers
fibthick = (self.hw - 2 * self.tf) / self.nw
area = fibthick * self.tw
yf = 0.5*self.hw - self.tf - 0.5*fibthick
for i in range(self.nw):
if self.mat[0] == 'Bilinear':
fiber = bilinearMaterial(self.mat)
elif self.mat[0] == 'BoucWen':
fiber = boucWenMaterial(self.mat)
else:
print("Error: Wrong material given to wide-flange section")
import sys
sys.exit()
self.fibers.append(fiber)
self.locations.append(yf)
self.areas.append(area)
yf = yf - fibthick
# -------------------------------------------------
# State determination
# -------------------------------------------------
def stateDetermination(self, us):
# Store total deformations to accommodate element13 & total strain materials
self.trialUsTotal0 = self.committedUsTotal0 + us[0, 1]
self.trialUsTotal1 = self.committedUsTotal1 + us[1, 1]
Ks = np.zeros((2, 2))
Fs = np.zeros(2)
for i in range(self.numfib):
# Fiber information
fiber = self.fibers[i]
area = self.areas[i]
loc = self.locations[i]
# Transformation
Tms = np.array([-loc, 1])
# Strain
epsilon = Tms.dot(us)
# Material state determination for fiber
sigma, E = fiber.stateDetermination(epsilon)
# Section integration
Ks = Ks + np.outer(Tms, Tms) * area * E
Fs = Fs + Tms * area * sigma
return Fs, Ks
# -------------------------------------------------
# DDM sensitivity analysis
# -------------------------------------------------
def stateDerivative(self, us, ddmParameter, ddmIndex, ddmIsHere, dKflag='none'):
dKs = np.zeros((2,2))
dFs = np.zeros(2)
dKsdus = np.zeros((2,2,2))
if ddmParameter in ['hw', 'bf', 'tf', 'tw'] and ddmIsHere:
# FLANGES
if ddmParameter == 'hw':
dFibthick = 0.0
dArea = 0.0
dyf = 0.5
elif ddmParameter == 'bf':
dFibthick = 0.0
fibthick = self.tf / self.nf
dArea = fibthick
dyf = 0.0
elif ddmParameter == 'tf':
dFibthick = 1.0 / self.nf
dArea = dFibthick * self.bf
dyf = -0.5 * dFibthick
elif ddmParameter == 'tw':
dFibthick = 0.0
dArea = 0.0
dyf = 0.0
for i in range(0, self.nf * 2, 2):
# TOP
fiber = self.fibers[i]
area = self.areas[i]
loc = self.locations[i]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
sigma, E = fiber.stateDetermination(epsilon)
dSigma, dE, dKmdum = fiber.stateDerivative(epsilon, ddmParameter, ddmIndex, False, dKflag)
dTms = np.array([-dyf, 0])
dEpsilon = dTms.dot(us)
dFs = dFs + dTms * area * sigma + Tms * dArea * sigma + Tms * area * E * dEpsilon[0] + Tms * area * dSigma
if dKflag == 'Initial':
E = fiber.initialStiffness()
dKs = dKs + np.outer(dTms, Tms) * area * E + np.outer(Tms, dTms) * area * E + np.outer(Tms, Tms) * dArea * E + np.outer(Tms, Tms) * area * (dKmdum * dEpsilon[0] + dE)
A, B, C = np.ix_(Tms, Tms, Tms)
dKsdus += np.multiply(dKmdum * area, A * B * C)
# BOTTOM
fiber = self.fibers[i+1]
area = self.areas[i+1]
loc = self.locations[i+1]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
sigma, E = fiber.stateDetermination(epsilon)
dSigma, dE, dKmdum = fiber.stateDerivative(epsilon, ddmParameter, ddmIndex, False, dKflag)
dTms = np.array([dyf, 0])
dEpsilon = dTms.dot(us)
dFs = dFs + dTms * area * sigma + Tms * dArea * sigma + Tms * area * E * dEpsilon[0] + Tms * area * dSigma
if dKflag == 'Initial':
E = fiber.initialStiffness()
dKs = dKs + np.outer(dTms, Tms) * area * E + np.outer(Tms, dTms) * area * E + np.outer(Tms, Tms) * dArea * E + np.outer(Tms, Tms) * area * (dKmdum * dEpsilon[0] + dE)
A, B, C = np.ix_(Tms, Tms, Tms)
dKsdus += np.multiply(dKmdum * area, A * B * C)
dyf = dyf - dFibthick
# WEB
if ddmParameter == 'hw':
dFibthick = 1.0 / self.nw
dArea = dFibthick * self.tw
dyf = 0.5 - 0.5 * dFibthick
elif ddmParameter == 'bf':
dFibthick = 0.0
dArea = 0.0
dyf = 0.0
elif ddmParameter == 'tf':
dFibthick = -2.0 / self.nw
dArea = dFibthick * self.tw
dyf = -1.0 - 0.5 * dFibthick
elif ddmParameter == 'tw':
fibthick = (self.hw - 2 * self.tf) / self.nw
dFibthick = 0.0
dArea = fibthick
dyf = 0.0
for i in range(2 * self.nf, 2 * self.nf + self.nw, 1):
fiber = self.fibers[i]
area = self.areas[i]
loc = self.locations[i]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
sigma, E = fiber.stateDetermination(epsilon)
dSigma, dE, dKmdum = fiber.stateDerivative(epsilon, ddmParameter, ddmIndex, False, dKflag)
dTms = np.array([-dyf, 0])
dEpsilon = dTms.dot(us)
dFs = dFs + dTms * area * sigma + Tms * dArea * sigma + Tms * area * E * dEpsilon[0] + Tms * area * dSigma
if dKflag == 'Initial':
E = fiber.initialStiffness()
dKs = dKs + np.outer(dTms, Tms) * area * E + np.outer(Tms, dTms) * area * E + np.outer(Tms, Tms) * dArea * E + np.outer(Tms, Tms) * area * (dKmdum * dEpsilon[0] + dE)
A, B, C = np.ix_(Tms, Tms, Tms)
dKsdus += np.multiply(dKmdum * area, A * B * C)
dyf = dyf - dFibthick
else:
for i in range(self.numfib):
fiber = self.fibers[i]
area = self.areas[i]
loc = self.locations[i]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
dSigma, dE, dKmdum = fiber.stateDerivative(epsilon, ddmParameter, ddmIndex, ddmIsHere, dKflag)
dKs = dKs + np.outer(Tms, Tms) * area * dE
dFs = dFs + Tms * area * dSigma
A, B, C = np.ix_(Tms, Tms, Tms)
dKsdus += np.multiply(dKmdum * area, A * B * C)
return dFs, dKs, dKsdus
# -------------------------------------------------
# Equilibrium iterations
# -------------------------------------------------
def equilibrium(self, Fs):
# Set iteration controls
maxiter = 100
tol = 1e-5
# Initialize trial displacement vector for iteration
usTrial = np.zeros(2)
# Iterate on section equilibrium
for j in range(maxiter):
# This method is called by element13, which works with incremental deformations.
# That means vsTrial is the second column of the U-matrix. Keep a history variable
# to give total deformation to materials that work with that.
self.trialUsTotal0 = self.committedUsTotal0 + usTrial[0]
self.trialUsTotal1 = self.committedUsTotal1 + usTrial[1]
FsTrial, Ks = self.stateDetermination(np.array([[self.trialUsTotal0, usTrial[0], 0.0],
[self.trialUsTotal1, usTrial[1], 0.0]]))
FsResidual = Fs - FsTrial
residualNorm = np.linalg.norm(FsResidual)
if residualNorm < tol:
break
usTrial = np.linalg.solve(Ks, FsResidual) + usTrial
if residualNorm >= tol:
print("ERROR: Equilibrium iterations in section failed with residual", residualNorm)
import sys
sys.exit()
return usTrial, Ks
# -------------------------------------------------
# Commit
# -------------------------------------------------
def commit(self):
# Commit total deformation to accommodate element13
self.committedUsTotal0 = self.trialUsTotal0
self.committedUsTotal1 = self.trialUsTotal1
# Commit material states
for i in range(self.numfib):
fiber = self.fibers[i]
fiber.commit()
# -------------------------------------------------
# Commit sensitivity history variables
# -------------------------------------------------
def commitSensitivity(self, us, ddmus, ddmParameter, ddmIndex, ddmIsHere):
if ddmParameter in ['hw', 'tf']:
# FLANGES
if ddmParameter == 'hw':
dFibthick = 0.0
dArea = 0.0
dyf = 0.5
elif ddmParameter == 'tf':
dFibthick = 1.0 / self.nf
dArea = dFibthick * self.bf
dyf = -0.5 * dFibthick
for i in range(0, self.nf * 2, 2):
# TOP
fiber = self.fibers[i]
loc = self.locations[i]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
dTms = np.array([-dyf, 0])
ddmepsilon = dTms.dot(us[:, 0]) + Tms.dot(ddmus)
fiber.commitSensitivity(epsilon, ddmepsilon, ddmParameter, ddmIndex, False)
# BOTTOM
fiber = self.fibers[i+1]
loc = self.locations[i+1]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
dTms = np.array([dyf, 0])
ddmepsilon = dTms.dot(us[:, 0]) + Tms.dot(ddmus)
fiber.commitSensitivity(epsilon, ddmepsilon, ddmParameter, ddmIndex, False)
dyf = dyf - dFibthick
# WEB
if ddmParameter == 'hw':
dFibthick = 1.0 / self.nw
dyf = 0.5 - 0.5 * dFibthick
elif ddmParameter == 'tf':
dFibthick = -2.0 / self.nw
dyf = -1.0 - 0.5 * dFibthick
for i in range(2 * self.nf, 2 * self.nf + self.nw, 1):
fiber = self.fibers[i]
loc = self.locations[i]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
dTms = np.array([-dyf, 0])
ddmepsilon = dTms.dot(us[:, 0]) + Tms.dot(ddmus)
fiber.commitSensitivity(epsilon, ddmepsilon, ddmParameter, ddmIndex, False)
dyf = dyf - dFibthick
else:
for i in range(self.numfib):
fiber = self.fibers[i]
loc = self.locations[i]
Tms = np.array([-loc, 1])
epsilon = Tms.dot(us)
ddmepsilon = Tms.dot(ddmus)
fiber.commitSensitivity(epsilon, ddmepsilon, ddmParameter, ddmIndex, ddmIsHere)
# -------------------------------------------------
# Print response
# -------------------------------------------------
def getResponse(self):
Fs = np.zeros(2)
# Loop over fibers
for i in range(self.numfib):
# Fiber information
fiber = self.fibers[i]
area = self.areas[i]
loc = self.locations[i]
# Transformation
Tms = np.array([-loc, 1])
# Material state determination for fiber
sigma = fiber.getResponse()
# Section integration
Fs = Fs + Tms * area * sigma
return Fs