程序_结构风荷载理论与Matlab计算_风荷载_风荷载计算_元胞自动机_matlab怎么做风_风

%exam_PS_model CFS.m %3D GEF随高度变化(位移、弯矩与剪力) clear,clc,% close all, %结构参数 H=180,%结构高 dz=9; %离散 z=(H:-dz:dz)';z(1)=z(1)-1e-5; b=5.6;%结构宽 m0=10686; %质量分布 det=0;% 质量系数 sz1=ones(size(2)); zH =z/H; %风场参数 alfa=0.15;%剖面指数 z0=0.1;%粗糙长度 Iu=1./log(z/z0);%各纵向湍流度 Iv =0.75*Iu;%横向湍流度 Lu=300*(z/200).^0.555 %纵向积分尺度 Lv=0.25*Lu;%横向积分尺度 tao =1; T=600; du =6.868; dv =9.434; Czu =7; Czv = 6.5; %模态参数 pha =2.15;%振型指数 n1=0.26;%一阶频率 ks1=0.005;%阻尼比 %气动参数 Cd=0.8;%阻力系数 Cs=0.28;%尾流系数 cxu=Cd;cyv=Cd;cys=Cs; %尾流参数 Ls=3; s = 0.19; Ba=sqrt(2)*Iu; % UH=20; %结构顶参考风速 Uz=UH*zH.^alfa;%风速分布 %解析式变量 id=1;ib=2;is=3;%响应类型e ndx =[id ib is]; Ju=2*Iu(1);Jv=Iv(1);Js=1; %Table 2 ma1=m0*H*(1/(2*pha +1)-det/(2*pha+2)); zHO = zH.pha;zH1=zH0.*zH; zH2=zH1.*zH;zH3 =zH2.*zH; mae(:,id)=zH0/(2*pi*nl)^2; mae(:,ib)=m0*H*H*((1+det*zH)/(pha+2).*(1-zH2)-... zH/(pha +1).*(1-zH1)-det/(pha +3)*(1-zH3)); mae(:,is)=m0*H*((1-zH1)/(pha +1)-det/... (pha+2)*(1-zH2)); %Table 3 zH6 =zH.^(alfa +1);zH7 =zH6.*zH; zH8=zH.^(2*alfa +1);zH9=zH8.*zH; Kp(:,id)=1/(ma1*(2*pi*n1)^2)/(alfa +pha +1)*zH0; Kp(:,ib)=H/(alfa +2)*(1-zH7)-z/(alfa +1).*(1-zH6); Kp(:,is)=1/(alfa+1)*(1-zH6); K(:,id)=1/(ma1*(2*pi*n1)^2)/(2*alfa +pha +1)*zH0; K(:,ib)=H/(2*alfa+2)*(1-zH9)-z/(2*alfa+1).*(1-zH8); K(:,is)=1/(2*alfa+1)*(1-zH8); % 式(44) zauv(:,id)=0.6*H*sz1; zauv(:,ib)=0.6*H+0.4*z; zauv(:,is)=0.6*H+0.4*z; zas=0.8*H*sz1; %式(45),(46),(47) kuv(:,id)=0.5*exp(-0.27*alfa*1)*exp(-0.27*pha)*sz1; ks(:,id)=0.5*exp(-0.27*alfa*2)*exp(-0.27*pha)*s21; kuv(:,ib)=0.5*(1-zH).*exp(-0.27*alfa*1*... exp(-2.3*zH))*exp(-0.27); ks (:,ib)=0.5*(1-zH).*exp(-0.27*alfa*2*... exp(-2.3*zH))*exp(-0.27); kuv(:,is)=0.5*(1-zH).*exp(-0.27*alfa*1*... exp(-3.5*zH)); ks(:,is)=0.5*(1-zH).*exp(-0.27*alfa*2*... exp(-3.5*zH)); %等效响应类型循环 cr='rbm';%plot line color for k=1:3 %u,v ips=ndx(k);%响应类型 ui=interpl(z,Uz,zauv(:,ips),'linear','extrap'); Lui=interpl(z,Lu,zauv(:,ips),'linear','extrap'); Lvi=interpl(z,Lv,zauv(:,ips),'linear','extrap'); %式(56)-(57) Tau(:,ips)=tao*du*ui/du./Lui; Tav(:,ips)=tao*du*ui/dv./Lvi; Lau(:,ips)=kuv (:,ips)*du*Czu*H/du./Lul; Lav(:,ips)=kuv(:,ips)*du*Czv*H/dv./vi; %式(53)-(55) Qu(:,ips)=1./(1 +0.56*(Tau(:,ips)).^0.74+... 0.3*(Lau(:,ips)).^0.63); Qv(:,ips)=1./(1 +0.56*(Tav(:,ips)).^0.74 +... 0.3*(Lav(:,ips)).^0.63); mvu(:,ips)=du*ui/du./Lui./sqrt(31.25*... (Tau(:,ips)).^1.44 +0.74*(Lau (:,ips)).^0.64 +... 5.41*(Tau(:,ips)).^0.93.*(Lau(:,ips)).^0.71); mvv(:,ips)=du*ui/dv./vl./sgrt(31.25*... (Tav(:,ips)).^1.44 +0.74*(Lav(:,ips)).^0.64+... 5.41*(Tav(:,ips)).^0.93.*(Lav(:,ips)).^0.71); ifk==1, ndu =n1*du*Lui/du./ui; ndv =n1*dv*Lvi/du./ui; wg=ndu.*Lau(:,id); Cf=1./wg-0.5*(1-exp(-2*wg))./wg.^2; Du =pi/4/ks1* du* ndu./... (1 +1.5* du* ndu).^ (5/3).* Cf; wg=ndv.*Lav(:,id); Cf =1. /wg-0.5* (1-exp (-2* wg))./wg.^2; Dv=pi/4/ks1* du* ndv./... (1 +1.5* du* ndv).^(5/3) .* Cf; end %s%式(59)-(64) Las(:,ips)=ks(:,ips)*H/Ls/b; if k==1, us=interp1(z,Uz,zas,'linear','extrap'); Bs =interpl(z,Ba,zas,'linear','extrap'); nas=nl*b/s./us; Fns=nas.^4./(nas.^4-3*nas.^2 +4); mvs =n1*Fns./nas; wg=Las(:,id); Cf=(1./wg-0.5*(1-exp(-2*wg))./wg.^2); Ds=sqrt(pi)*nas/4/ks1./Bs.*... exp(-((1-nas)./Bs).^2).*Cf; end wg=Las(:,ips); Cf =1./wg-0.5*(1-exp(-2*wg))./wg.^2; Qs(:,ips)=Cf.*Fns; %式(33) Xxu(:,ips)=Ju*Kp(:,ips)./R(:,ips)*cxu/cxu; Xyv(:,ips)=Jv*Kp(:,ips)./K(:,ips)*cyv/cxu; Xys(:,ips)=Js*K(:,ips)./K(:,ips)*cys/cxu; %式(34) Fxu(:,ips)=Ju*Kp(:,id)./K(:,ips).*... mae(:,ips)./mae(:,id)*cxu/cxu; Fyv(:,ips)=Jv*Kp(:,id)./K(s,ips).*... mae(:,ips)./mae(:,id)*cyv/cxu; Fys(:,ips)-Js*R(:,id)./K(:,ips).*... mae(:,ips)./mae(:,id)*cys/cxu %式(30) Qx(:,ips)=xxu(:,ips).^2.*Qu(:,ips); Qy(:,ips)=xyv(:,ips).^2.*Qvt(:,ips)+... xys(:,ips).^2.*Qs(:,ips); %式(31) Dx(:,ips)=Fxu(:,ips).^2.*Du(:,id); Dy(:,ips)=Fyv(:,ips).^2.*Dv(:,id)+... Fys(:,ips).^2.*Ds(:,id); %式(32) mvx(:,ips)=Xxu(:,ips).^2.*Qu(:,ips).*... mvu(:,ips)./Qx(:,ips); mvy(:,ips)=(xyv(:,ips).^2.*Qv(:,ips).*... mvv(:,ips)+Xys(:,ips).^2.*... Qs(:,ips).*mvs)./Qy(:,ips); %式(25) vvx(:,ips)=sqrt((mvx(:,ips).*Qx(:,ips)+... n1^2*Dx(:,ips))./(Qx(:,ips)+Dx(:,ips))); vvy(:,ips)=sqrt((mvy(:,ips).*gy(:,ips)+... n1^2*Dy(:,ips))./(Qy(:,ips)+Dy(:,ips))); %式(5) tmpx =sqrt(2*log(T*vvx(:,ips))); grx(:,ips)=tmpx +0.5772./tmpx; tmpy=sqrt(2*log(2*T*wvy(:,ips))); gry(:,ips)=tmpy +0.5772./tmpy; % %3D GEF Gax(:,ips)=1 +grx(:,ips).*... sqrt(Qx(:,ips)+Dx(:,ips)); Gay(:,ips)=gry(:,ips).*... sqrt(Qy(:,ips)+Dy(:,ips)); % %show figure(1),hold on plot(Gax(:,ips),zH,[cr(k),'-']) hold off figure(2),hold on, plot(Gay(:,ipa),zH,[cr(k),'-']) hold off % %end return