Ruben Jerves 9/8/2023
Homework 2
CEE 5470
Problem 1)
Ground motion
R = 0.05; Ti = 0.1; Tf = 3.0; dT = 0.01;
Spectro = ResponseSpectra(R,Ti,Tf,dT);
Pseudo-acceleration response spectrum
x = Spectro(:,1); y = Spectro(:,4);
x_label = 'Natural perion, Tn (s) ';
y_label = 'Sa (g)';
title_label = 'Pseudo-acceleration response spectrum';
Graph(x,y,x_label,y_label,title_label)
Page 1 of 11
Ruben Jerves 9/8/2023
Displacement response spectrum
x = Spectro(:,1); y = Spectro(:,2);
x_label = 'Natural perion, Tn (s) ';
y_label = 'Sd (in)';
title_label = 'Displacement response spectrum';
Graph(x,y,x_label,y_label,title_label)
Page 2 of 11
Ruben Jerves 9/8/2023
Problem 2)
% Units
inch = 1;
kips = 1;
s = 1;
ft = 12*inch;
g = 386.1*inch/(s^2);
% Data
E = 29000 *kips/(inch^2);
I = 28.1 * inch^4;
h = 12* ft;
d = 6.625 *inch;
Calculations
M = 30*kips/g % Mass [kip-s^2/in]
M = 0.0777
Ksys = 6 * (3*E*I)/(h^3) % System stiffness [kip/in]
Ksys = 4.9124
wn = (Ksys/M).^0.5 % Natural frequency [rad/s]
wn = 7.9512
Tn = 2*pi()/wn % Natural perio [s]
Tn = 0.7902
delta_max = 2.8 *inch % Peak displacement (from spectrum)
delta_max = 2.8000
V_max = M*(delta_max*wn^2) % Max. shear in the columns (can be also obtained
the spectum of pseudo-acceleration) [kips]
V_max = 13.7546
M_max = V_max*h % Max. moment [kip-in]
M_max = 1.9807e+03
Sigma_max = (M_max*(0.5*d))/I %Max bending stress at base of columns [ksi]
Sigma_max = 233.4852
Page 3 of 11
Ruben Jerves 9/8/2023
Problem 3 )
h = 14 *ft; % story height
N_stories = 5; % Number of stories
H_building = N_stories*h % Building height [in]
H_building = 840
Ta = 0.016*(H_building/12)^0.9 % Approximated period, based on ASCE 7-22
(12.8.8)
Ta = 0.7323
global Base_shear Roof_disp
Base_shear = [0,300,600,780,880,900,880,852,800,720,600];
drift = (0:0.5:5)*(1/100);
Page 4 of 11
Ruben Jerves 9/8/2023
Roof_disp = drift*H_building;
Graph(Roof_disp,Base_shear,'Roof displacement, \Delta (in)','Base shear, V
(kips)','Building pushover')
a.
f(812)
ans = 38.1920
c = fminbnd(@f,0,1000) % this would be the optimized (minumun); in other words
Vy
c = 810.7535
val = f(c) % For instance, here it is seen how f(c) was minimized
val = 5.3470e-04
Ko = (Base_shear(2) - Base_shear(1))/(Roof_disp(2) - Roof_disp(1)); % Initial
stiffness [kip/in]
x_inter = c/Ko % This would be the delta at yielding
x_inter = 11.3505
x1 = 0:0.1:x_inter;
y1 = Ko*x1;
x2 = x_inter:0.1:Roof_disp(end);
y2 = c*ones(length(x2),1);
Page 5 of 11
Ruben Jerves 9/8/2023
figure()
hold on
grid on
plot(x1,y1,'LineWidth',1.5,'Color','blue');
plot(x2,y2,'LineWidth',1.5,'Color','blue');
plot(Roof_disp,Base_shear,'LineWidth',1.5,'Color','red')
title('Equal energy elasto-plastic idealization','FontSize',12);
ylabel('Base shear, V (kips)','FontSize',12);
xlabel('Roof displacement (in)','FontSize',12);
axes_handle=get(gcf,'CurrentAxes');
set(axes_handle,'FontSize',12,'FontName','Arial');
set(gcf,'Color',[1,1,1])
grid on
Another way to obtain these values, Vy and delta y, can be as fallow:
I assume that the max displacement, would be the one that corresponds to 80% of the peak load,
thus:
delta_m = 37.8 * inch
delta_m = 37.8000
V_yield = 840 * kips % Base shear at yielding (based on graph)
V_yield = 840
Page 6 of 11
Ruben Jerves 9/8/2023
delta_yield = 12 * inch % displacement at yielding (based on graph)
delta_yield = 12
For instance, the base shear and displacement at yielding are obtained from the graph, by
approximating the areas, . This procedure was done graphically.
b.
The ductility of the system would be as follows:
delta_y = x_inter
Page 7 of 11
Ruben Jerves 9/8/2023
delta_y = 11.3505
delta_u = Roof_disp(end)
delta_u = 42
mu_delta = delta_u/delta_y % This would be our mu
mu_delta = 3.7003
Now, the ductility factor R_mu is calculated as fallow:
Ta
Ta = 0.7323
R_mu_1 = (2*mu_delta-1)^0.5
R_mu_1 = 2.5299
R_mu_2 = mu_delta
R_mu_2 = 3.7003
Our value of R_mu would fall between R)mu_1 and R_mu_1; thus, we interpolate to find R_mu.
R_mu = interp1([0.5, 1.0],[R_mu_1,R_mu_2],Ta) % This would be our R_mu
R_mu = 3.0737
c.
The overstrength factor is calculated as fallow:
V_y = c
V_y = 810.7535
V_design = c/3
V_design = 270.2512
R_s= max(Base_shear) / V_design
R_s = 3.3302
d.
The appropriate R factor is as follows:
R_r = 1.0;
Page 8 of 11
Ruben Jerves 9/8/2023
R = R_s*R_mu*R_r
R = 10.2363
Problem 4)
When using the overstrength factor as the sole criterion for designing structural elements, it may
result in inadequate performance under seismic loads, as it might not account for the higher
demands that a structural member could experience. For example, in situations where two or more
elements from the lateral force resisting system converge at a single point, if one of these elements
buckles under compression, it can generate additional moments around that node coming from the
elements in tension. This scenario is particularly relevant in the context of concentrically braced
frames. Generally speaking, this can be usually the case when buckling can occur.
Functions
function Graph(x,y,x_label,y_label,title_label)
figure()
hold on
grid on
plot(x,y,'LineWidth',1.5);
title(title_label,'FontSize',12);
ylabel(y_label,'FontSize',12);
xlabel(x_label,'FontSize',12);
axes_handle=get(gcf,'CurrentAxes');
set(axes_handle,'FontSize',12,'FontName','Arial');
set(gcf,'Color',[1,1,1])
grid on
end
function Spectro = ResponseSpectra(R,Ti,Tf,dT)
close all
nTn = ceil(((Tf-Ti)/dT));
Tn = linspace(Ti,Tf,nTn);
dt = 0.01;
u0 = 0;
v0 = 0;
gamma=0.5; % Newmark Coefficient
beta= 1/6; % Newmark Coefficient (linear acceleration method)
g=386.1; % [in/s2]
accel=load('Homework_2_NR94stn_acc.txt');
Page 9 of 11
Ruben Jerves 9/8/2023
t_end = dt*length(accel)-dt;
time_accel = 0:dt:t_end;
vg_dt=cumtrapz(accel,1)*dt;
ug_dt=cumtrapz(vg_dt,1)*dt;
for i = 1:nTn
wn = (2*pi())/Tn(i);
k = 20;
m = k/(wn^2);
Peff = -m*accel*g;
c = 2*m*R*wn;
[ur,~,ar,ut,~,at] =
NewmarksMethod(beta,gamma,c,m,k,dt,u0,v0,Peff,accel,vg_dt,ug_dt);
urmax = max([max(ur) abs(min(ur))]); utmax = max([max(ut)
abs(min(ut))]);
armax = max([max(ar) abs(min(ar))]); atmax = max([max(at)
abs(min(at))]);
Spectro(i,1) = Tn(i);
Spectro(i,2) = urmax; %Relative displacement (u or
Sd) --> used for structure forces
Spectro(i,3) = utmax; %Total displacement (u + ug)
Spectro(i,4) = (urmax*wn^2)/g; %Pseudo aceleracion
(Sd*wn^2/g)
end
Graph(time_accel',accel,'Time (s)','Ground acceleration (g)','Ground
motion')
end
function [ur,vr,ar,ut,vt,at] =
NewmarksMethod(beta,gamma,c,m,k,dt,u0,v0,p,ag_dt,vg_dt,ug_dt)
ur = []; vr = []; ar = [];
%Initial calculations
i = 1;
ur(i,1) = u0;
vr(i,1) = v0;
ar(i,1) = (p(i,1)-c*vr(i,1)-k*ur(i,1))/m;
a1 = (1/(beta*dt^2))*m + (gamma/(beta*dt))*c;
a2 = (1/(beta*dt))*m + ((gamma/beta)-1)*c;
Page 10 of 11
Ruben Jerves 9/8/2023
a3 = ((1/(2*beta))-1)*m + dt*((gamma/(2*beta))-1)*c;
k_gorro =k + a1;
%Calculations for each time step
n = length(p);
for i = 1:n-1
p_gorro(i+1,1) = p(i+1,1) + a1*ur(i,1) + a2*vr(i,1) + a3*ar(i,1);
ur(i+1,1) = p_gorro(i+1)/k_gorro;
vr(i+1,1) = (gamma/(beta*dt))*(ur(i+1,1)-ur(i,1)) + (1-
(gamma/beta))*vr(i,1) + dt*(1-(gamma/(2*beta)))*ar(i,1);
ar(i+1,1) = (1/(beta*dt^2))*(ur(i+1,1)-ur(i,1))-(1/(beta*dt))*vr(i,1)-
((1/(2*beta))-1)*ar(i,1);
end
at = ar + ag_dt;
vt = vr + vg_dt;
ut = ur + ug_dt;
end
% function to minimize (equal energy)
function val = f(c)
global Base_shear Roof_disp
Ko = (Base_shear(2) - Base_shear(1))/(Roof_disp(2) - Roof_disp(1)); %
Initial stiffness [kip/in]
x_inter = c/Ko;
A_tri = 0.5 * x_inter * c;
A_rec = (Roof_disp(end) - x_inter) * c;
A_total = A_tri + A_rec;
Area_under_curve = trapz(Roof_disp,Base_shear); % Area under the force vs
displacement curve (pushover) [kip-in]
val = abs(Area_under_curve - A_total);
end
Page 11 of 11
