matlab在建筑抗震工程中的应用 [大家看看,上个类似的帖子有误]
某三层钢筋混凝土结构，地震波采用200gal El Centro波，采样周期为0.02秒。用振型分解法求解结构地震反应的MATLAB程序如下:xs=2*0.287; % 地震波数据 dzhbo=load('ELCENTRO.txt', 'r'); ag=dzhbo*0.01*xs;dt=0.02; ndzh=400; cn=3; % 结构参数 m0=*1e+3; k0=*1e+5; l=diag(ones(cn));m=diag(m0); =matrixju(k0,cn); =eig(ik,m); % 结构动力特性的求解 d=diag(sqrt(d)); for i=1:cn =min(d);xgd(:,i)=x(:,j);d(j)=max(d)+1; end w=d1; x=xgd; a1=2*w(1)*w(2)*(0.05*w(2)-0.07*w(1))/(w(2)^2-w(1)^2);a2=2*(0.07*w(2)-0.05*w(1))/(w(2)^2-w(1)^2);for j=1:cn x(:,j)=x(:,j)/x(cn,j); znb0(j)=(a1+a2*w(j)^2)/2/w(j); zhcan(j)=(x(:,j))'*m*l/((x(:,j))'*m*x(:,j)); % 求解振型参与系数 =zxzj(znb0(j),w(j),ag); % 求解 file:///C:/Users/SHISHI~1/AppData/Local/Temp/msohtmlclip1/01/clip_image002.gif和file:///C:/Users/SHISHI~1/AppData/Local/Temp/msohtmlclip1/01/clip_image004.gifend for i=1:cn % 求解结构各层的地震反应 disp1=0; accel1=0; for j=1:cn disp0=zhcan(j)*dlt(j,:)*x(i,j);accel0=zhcan(j)*dltacceler(j,:)*x(i,j);disp1=disp1+disp0;accel1=accel1+accel0;end disp(i,:)=disp1;accel(i,:)=accel1;end t=0:dt:ndzh*dt; subplot(2,2,1);plot(t,disp(3,:)*1e+3,'k-');subplot(2,2,2);plot(t,accel(3,:),'k-');% This sub-program is solvingthe dynamic responses of single degree system function=zxzj(znb,w,dag) dt=0.02; n=400; x(1)=0; dx(1)=0; ddx(1)=0; s=1+znb*dt*w+w^2*dt^2/6; %中间参数 sfori=1:n a(i)=x(i)+dx(i)*dt+ddx(i)*dt^2/3;b(i)=dx(i)+ddx(i)*dt/2;ddx(i+1)=-1*(dag(i+1)+2*znb*w*b(i)+w^2*a(i))/s; % 加速度 dx(i+1)=b(i)+ddx(i+1)*dt/2; % 速度 x(i+1)=a(i)+ddx(i+1)*dt^2/6; % 位移 endbx=x;acceler=ddx;% This sub-program is formatrix aggregation of system function=matrixju(korc,cn) kcju=zeros(cn); for i=1:cn-1kcju(i,i)=korc(i)+korc(i+1); kcju(i,i+1)=-korc(i+1); kcju(i+1,i)=-korc(i+1); end kcju(cn,cn)=korc(cn)
运行之后 出现下面提示 ??? Error using ==> loadUnable to read file ELCENTRO.txt: No such file or directory.还有蓝色这个函数matrixju不知道什么意思？请各位大侠指教。
不是都说了
ELCENTRO.txt: No such file or directory
建议help下load 谢谢！晓得了，是我没要把文件放在工作目录里面，随便放在桌面山了 本人刚刚开始学习，以后请大侠多多指教 本帖最后由 牛小贱 于 2014-3-30 20:00 编辑
我把最后两个子程序存为matrixju.m和zxzj.m两个m文件后，程序如下，但是所得结果与已知结果不同，郁闷啊！！请大侠再指教！！不胜感激！！
xs=2*0.287; % 地震波数据
dzhbo=load('ELCENTRO.txt');
ag=dzhbo*0.01*xs;
dt=0.02;
ndzh=400;
cn=3; % 结构参数
m0=*1e+3;
k0=*1e+5;
l=diag(ones(cn));
m=diag(m0);
=matrixju(k0,cn);
=eig(ik,m); % 结构动力特性的求解
d=diag(sqrt(d));
for i=1:cn
=min(d);
xgd(:,i)=x(:,j);
d(j)=max(d)+1;
end
w=d1;
x=xgd;
a1=2*w(1)*w(2)*(0.05*w(2)-0.07*w(1))/(w(2)^2-w(1)^2);
a2=2*(0.07*w(2)-0.05*w(1))/(w(2)^2-w(1)^2);
for j=1:cn
x(:,j)=x(:,j)/x(cn,j);
znb0(j)=(a1+a2*w(j)^2)/2/w(j);
zhcan(j)=(x(:,j))'*m*l/((x(:,j))'*m*x(:,j)); % 求解振型参与系数
=zxzj(znb0(j),w(j),ag); % 求解 ∆_j和∆ ̈_j
end
for i=1:cn % 求解结构各层的地震反应
disp1=0;
accel1=0;
for j=1:cn
disp0=zhcan(j)*dlt(j,:)*x(i,j);
accel0=zhcan(j)*dltacceler(j,:)*x(i,j);
disp1=disp1+disp0;
accel1=accel1+accel0;
end
disp(i,:)=disp1;
accel(i,:)=accel1;
end
t=0:dt:ndzh*dt;
subplot(2,2,1)
plot(t,disp(3,:)*1e+3,'b-')
subplot(2,2,2)
页:
[1]