Hi, when I executed your code and my understanding is that the error is due to it is trying to access an data which is not computed yet.
for i = 2:1001
p_hat(i) = P(i)+(a1*u(i-1))+(a2*du(i-1))+(a3*ddu(i-1));
u(i) = p_hat(i)/k_hat;
du(i) = ((G/(B*delta_t))*(u(i)-u(i-1)))+((1-(G/B))*(du(i-1)))+(delta_t*(1-(G/(2*B)))*ddu(i-1));
ddu(i) = ((1/(B*(delta_t^2)))*(u(i)-u(i-1)))-((1/(B*delta_t))*du(i-1))-(((1/(2*B))-1)*ddu(i-1));
end
Here in this loop for the computation of “P_hat”, “u ” is needed which is evaluated in next expression.
However, you can initialize the values to get the plot.
figure
hold on
for m = 1000;
Z = 0.05;
t = 0:0.01:10;
wn = 0.5/Z;
k =(wn^2)*m;
G = 0.5;
B = (1/6);
wd = wn*(sqrt(1-(Z^2)));
delta_t = 0.01;
c = 2*Z*sqrt(k*m);
P = (100*sin(0.5*t))-(150*cos(0.5*t));
a1 = (m/(B*delta_t^2))+(G*c/(B*delta_t));
a2 = (m/(B*delta_t))+((G/B)-1)*c;
a3 = ((1/(2*B))-1)*m+delta_t*((G/(2*B))-1)*c;
k_hat = k + a1;
u=zeros(1,1000);
du=zeros(1,1000);
ddu=zeros(1,1000);
for i = 2:1001
p_hat(i) = P(i)+(a1*u(i-1))+(a2*du(i-1))+(a3*ddu(i-1));
u(i) = p_hat(i)/k_hat;
du(i) = ((G/(B*delta_t))*(u(i)-u(i-1)))+((1-(G/B))*(du(i-1)))+(delta_t*(1-(G/(2*B)))*ddu(i-1));
ddu(i) = ((1/(B*(delta_t^2)))*(u(i)-u(i-1)))-((1/(B*delta_t))*du(i-1))-(((1/(2*B))-1)*ddu(i-1));
end
plot(u)
xlabel('Time (sec)')
ylabel('Displacement (m)')
title('Time vs Displacement (Newmark Method)')
grid on
end
