jtao=[];jnr1=[];jnr2=[];
jnr71=[];jnr72=[];
jnr61=[];jnr62=[];
jnr51=[];jnr52=[];
jnr41=[];jnr42=[];
jnr31=[];jnr32=[];
for uii=1:5
ui=uii
ndzb=80000;ui=(ui-1)*ndzb+1;
sunweizhi=11;sunshangchengdu=0.002*(uii-1);
dt=0.005;
deta=0.50;arfa=0.25;
load ab_data.mat
% li=a(ui:ui+ndzb-1);clear acce;
li=a;
cn=21;
l=1.05/cn;E=2.06e11;b=0.04;h=0.005;rou=7850;b1=b-sunshangchengdu;
Iz=b*h*h*h/12;area=b*h;Iz1=b1*h*h*h/12;area1=b1*h;
k1=12/l/l;k2=6/l;k3=-k1;k4=-k2;
k0=[k1 k2 k3 k2;k2 4 k4 2;k3 k4 k1 k4;k2 2 k4 4];
ksun0=k0*E*Iz1/l;k0=k0*E*Iz/l;
m1=156;m2=22*l;m3=54;m4=-13*l;
m5=4*l*l;m6=13*l;m7=-3*l*l;m8=-22*l;
m0=[m1 m2 m3 m4;m2 m5 m6 m7;m3 m6 m1 m8;m4 m7 m8 m5];
msun0=rou*area1*l/420*m0;
m0=rou*area*l/420*m0;
K=zeros(cn*2+2);M=K;C=K;
for i=1:cn
j=2*i-1;
if i==sunweizhi
K(j:j+3,j:j+3)=K(j:j+3,j:j+3)+ksun0;
M(j:j+3,j:j+3)=M(j:j+3,j:j+3)+msun0;
else
K(j:j+3,j:j+3)=K(j:j+3,j:j+3)+k0;
M(j:j+3,j:j+3)=M(j:j+3,j:j+3)+m0;
end
end
K(1:2,:)=[];K(:,1:2)=[];
M(1:2,:)=[];M(:,1:2)=[];
a0=1/(arfa*dt*dt);a1=deta/(arfa*dt);a2=1/(arfa*dt);a3=0.5/arfa-1;a4=deta/arfa-1;a5=dt/2*(deta/arfa-2);
a6=dt*(1-deta);a7=deta*dt;
ch=zeros(cn*2,1);
wy1=ch;sd1=ch;jsd1=ch;
ch=zeros(7,ndzb);
wyt=ch;
unit=ones(cn*2,1);
for i=1:cn
unit(i*2)=0;
end
[x,d1]=eig(K,M);d=sqrt(d1);w=sort(diag(d));
a=2*w(1)*w(2)*(0.1*w(2)-0.1*w(1))/(w(2)^2-w(1)^2);
b=2*(0.1*w(2)-0.1*w(1))/(w(2)^2-w(1)^2);
C=a*M+b*K;
ksin=K+a1*C+a0*M;
kni=inv(ksin);
for i=1:ndzb
pxin=-M*unit*li(i)+M*(a0*wy1+a2*sd1+a3*jsd1)+C*(a1*wy1+a4*sd1+a5*jsd1);
xdt=kni*pxin;
jsd=a0*(xdt-wy1)-a2*sd1-a3*jsd1;jsd0=jsd+li(i);
sd=sd1+a7*jsd+a6*jsd1;
wy1=xdt;sd1=sd;jsd1=jsd;
for jj=1:7
wyt(jj,i)=wy1(jj*6-1);
end
end
s=wyt;
save li.txt li -ASCII -TABS
save js.txt s -ASCII -TABS
kik=0;
kl=10;nn=1000;kkk=10000;xwai=5;
load js.txt;s=js';load li.txt;
mn=7;
mm=nn/4;
dk=nn;%数据增量
haj=zeros(nn,kl,mn);hjp=zeros(mn,nn);
for i1=1:kl
k=kkk+(i1-1)*dk;fj=lizhen(li,nn,k);aa=fj';aa=inv(aa); i1
for zz=1:mn
a=s(k:k+nn-1,zz);x=aa*a;haj(:,i1,zz)=x';hjp(zz,:)=hjp(zz,:)+x';
% kik=kik+1;figure(kik);plot(x);
end
end
hjp=hjp/kl;
wch2=zeros(mn,kl);
wch=zeros(mn,kl);wc=ones(1,mn)*100;wc1=wc;wch1=zeros(mn,1);
for zz=1:mn
for k1=1:kl
for i=2:nn
wch(zz,k1)=wch(zz,k1)+abs(haj(i,k1,zz)-haj(i-1,k1,zz));
end
wc(zz)=min(wc(zz),wch(zz,k1));
x=haj(:,k1,zz);
wch2(zz,k1)=sum(abs(x-hjp(zz,:)'));
wch1(zz)=wch1(zz)+wch2(zz,k1);
end
end
wch1=wch1/kl;
ni=zeros(1,mn);lo=zeros(mn,kl);wchh=wch;wc1=wc;
for zz=1:mn
for k1=1:xwai
[xx,yy]=min(wchh(zz,:));
if (xx<wc1(zz)*11) %(wch2(zz,yy)<wch1(zz)*0.7)
wchh(zz,yy)=10000;
lo(zz,yy)=yy;
ni(zz)=ni(zz)+1;
wc(zz)=xx;
end
end
end
hp=zeros(mn,nn);
for zz=1:mn
for i1=1:kl
if lo(zz,i1)>0
hp(zz,:)=hp(zz,:)+haj(:,i1,zz)';
end
end
hp(zz,:)=hp(zz,:)/ni(zz);
end
save jhp.txt hp -ASCII -TABS
sk=zeros(mm,mn,mn);sk1=zeros(mn,mn);ssk=sk;
for ii=1:mn
for jj=1:mn
sk1(ii,jj)=huxg1(hp(ii,:),hp(jj,:),nn,0);
end
end
load mtao.txt;load mnr1.txt;load mnr2.txt;
tao=zeros(mn,mn);nr1=tao;nr2=tao;
for ii=1:mn
for jj=1:mn
nr1(ii,jj)=huxg1(hp(ii,:),hp(jj,:),nn,mtao(ii,jj))/sk1(jj,jj);
nr2(ii,jj)=huxg1(hp(ii,:),hp(jj,:),nn,mtao(ii,jj))/sk1(ii,jj);
end
end
jnr1=[jnr1;nr1];jnr2=[jnr2;nr2];
jtao=[jtao;tao];
save jnr1.txt jnr1 -ASCII -TABS
save jtao.txt jtao -ASCII -TABS
save jnr2.txt jnr2 -ASCII -TABS
% if uii==1
% save mnr1.txt nr1 -ASCII -TABS
% save mtao.txt tao -ASCII -TABS
% save mnr2.txt nr2 -ASCII -TABS
% end
load mnr1.txt;
load mnr2.txt;
%差值
nr31=nr1-mnr1;jnr31=[jnr31;nr31];
nr32=nr2-mnr2;jnr32=[jnr32;nr32];
save jnr31.txt jnr31 -ASCII -TABS
save jnr32.txt jnr32 -ASCII -TABS
%变化率
nr41=nr31./mnr1*100;jnr41=[jnr41;nr41];
nr42=nr32./mnr2*100;jnr42=[jnr42;nr42];
save jnr41.txt jnr41 -ASCII -TABS
save jnr42.txt jnr42 -ASCII -TABS
%绝对值
nr51=abs(nr1);jnr51=[jnr51;nr51];
nr52=abs(nr2);jnr52=[jnr52;nr52];
save jnr51.txt jnr51 -ASCII -TABS
save jnr52.txt jnr52 -ASCII -TABS
%绝对值差值
nr61=nr51-abs(mnr1);jnr61=[jnr61;nr61];
nr62=nr52-abs(mnr2);jnr62=[jnr62;nr62];
save jnr61.txt jnr61 -ASCII -TABS
save jnr62.txt jnr62 -ASCII -TABS
%绝对值变化率
nr71=nr61./abs(mnr1)*100;jnr71=[jnr71;nr71];
nr72=nr62./abs(mnr2)*100;jnr72=[jnr72;nr72];
save jnr71.txt jnr71 -ASCII -TABS
save jnr72.txt jnr72 -ASCII -TABS
%
%
%
% fg=zeros(mn,mn);
% tao=10;
% flag=min(min(fg));
% while flag==0
%
% for jj=1:mn
% slin=sk(:,ii,jj);ji=slin(tao)*slin(tao+1);
% if (ji<=0)&(fg(ii,jj)==0)
% if abs(slin(tao))<=abs(slin(tao+1))
% tao0=tao;
% else
% tao0=tao+1;
% end
% tao1=tao+slin(tao)/(slin(tao)-slin(tao+1));m3(ii,jj)=tao1-1;
% m4(ii,jj)=(slin(tao)-slin(tao+1))*(tao+1-tao1)+slin(tao+1);
% m1(ii,jj)=slin(tao0);m2(ii,jj)=tao0-1;fg(ii,jj)=1;
% end
% end
% end
% flag=min(min(fg))
% tao=tao+1;
% end
end
没有合适的资源?快使用搜索试试~ 我知道了~
温馨提示
【达摩老生出品,必属精品,亲测校正,质量保证】 资源名:matlab_分析机械力学中的工程振动问题 资源类型:matlab项目全套源码 源码说明: 全部项目源码都是经过测试校正后百分百成功运行的,如果您下载后不能运行可联系我进行指导或者更换。 适合人群:新手及有一定经验的开发人员
资源推荐
资源详情
资源评论
收起资源包目录
matlab_分析机械力学中的工程振动问题.rar (1个子文件)
good0412.m 6KB
共 1 条
- 1
资源评论
- m0_697162112022-04-25用户下载后在一定时间内未进行评价,系统默认好评。
- 工控老马2022-05-26用户下载后在一定时间内未进行评价,系统默认好评。
阿里matlab建模师
- 粉丝: 3334
- 资源: 2784
下载权益
C知道特权
VIP文章
课程特权
开通VIP
上传资源 快速赚钱
- 我的内容管理 展开
- 我的资源 快来上传第一个资源
- 我的收益 登录查看自己的收益
- 我的积分 登录查看自己的积分
- 我的C币 登录后查看C币余额
- 我的收藏
- 我的下载
- 下载帮助
安全验证
文档复制为VIP权益,开通VIP直接复制
信息提交成功