脉动风时程matlab程序

B类风场与台风风场下输电塔的风振响应和风振系数楼文娟;夏亮;蒋莹;金晓华;王振华【摘要】为研究常规B类风场与台风风场下输电塔的风振响应差异,以沿海地区某四回路角钢输电塔为原型建立了有限元模型,采用谐波叠加法生成两类风场下的风速时程,并在时域内进行了输电塔风振响应和风振系数的数值分析.结果表明:台风风场的高湍流特性导致其作用时各测点的顺风向风振响应均大于B类风场下的对应值.两类风场下,输电塔的风振系数比值约为1.25.因此,台风多发地区的输电塔设计必须考虑台风高湍流引起的动力风荷载增大效应.此外还进行了气弹模型风洞试验,以研究不同风速下的气动响应和风振系数,并将试验结果与理论计算进行了分析比较,验证了数值分析的适用性.%A numerical analysis on the wind-induced response of a four-circuit angle-steel transmission tower under conventional terrain B wind field and typhoon wind field was performed. A FEM model was established, and the dynamic response was calculated under a fluctuating wind field simulated by using harmonic wave superimposing method. Based on data of numerical analysis, wind-induced responses under each wind field were discussed. Essential conclusions are as follows; high turbulence and strong variability of typhoon wind field have great influence on the response of transmission towers. The RMS of acceleration under typhoon wind field is larger than that under terrain B wind field. Under the two types of wind fields, the average ratio of wind load factor is about 1.25. Therefore, the design of transmission towers in typhoon-prone areas should take the fluctuating wind load magnification effect into consideration. Furthermore, the wind tunnel test on anaeroelastic model of the transmission tower was performed to study its wind-induced responses under different velocity. The test results were compared with theoretical values and the accuracy of the numerical analysis was verified.【期刊名称】《振动与冲击》【年(卷),期】2013(032)006【总页数】5页(P13-17)【关键词】输电塔;数值分析;风振响应;风振系数;台风风场【作者】楼文娟;夏亮;蒋莹;金晓华;王振华【作者单位】广东省电力设计研究院,广州510663【正文语种】中文【中图分类】TU973.32我国东南沿海为台风多发地区，台风风场的高湍流度、强离散性和强变异性等特征将产生与良态风作用下不同的复杂风振效应，而现行规范尚未涉及台风作用下输电塔风荷载的具体规定。

function p=pulsegen(fs,T,edge,type,f,opt);%p=pulsegen(fs,T,edge,type,f,opt);%a signal generation program%fs is the sampling frequency%T is the total signal length%edge is a decay parameter for some waveforms% it is used in 'gaussian', 'monocycle', 'biexponential', 'mexican hat', 'sinc', 'double sinc', 'sinc squared'% and windowed sweep% it is mostly a parameter to describe how much the edge of the pulse is decayed.%type is the type of the waveform desired% allowable types are 'gaussian', 'square', 'triangle', 'monocycle',% 'biexponential', 'mexican hat', 'sinc', 'double sinc', 'sinc squared','sweep', and 'windowed sweep'%f is the modulation frequency if left out it is assumed 0.%opt is an optional argument for pulse waveforms requiring a lower and higher frequency% it is used in 'double sinc' ,'sweep' and 'windowed sweep' for the low and high frequency.% the pulses are always normalized to a peak amplitude of 1if nargin<4 %test for optional argumentserror('not enough input arguments');elseif nargin==4f=0;opt=[16*edge/(5*T),64*edge/(5*T)];elseif nargin==5opt=[16*edge/(5*T),64*edge/(5*T)];endif (edge==0)edge==1;endt=-T/2:1/fs:T/2;sig=(T/8/edge)^2;switch typecase {'guassian'} %generate a guassian pulsey=exp(-(t).^2/sig);case {'square'} %generate a square pulsey=ones(size(t));case {'triangle'} %generate a triangle pulsey=(t+T/2).*(t<0)-(t-T/2).*(t>=0);case {'monocycle'} %generate a gaussian monocycley=2*t./sig.*exp(-(t).^2/sig);case {'biexponential'} %generate a biexponential pulsey=exp(-abs(t)*8*edge/T);case {'mexican hat'} %generate a gaussian second deriviativez=t./sqrt(0.75*sig);y=sqrt(1/2*pi).*(1-z.^2).*exp(-z.^2/2);case {'sinc'} %generate a sinc functiony=sinc(2*pi*edge*16.*t/(5*T));case {'double sinc'} %generate a bandlimited function from two sinc functions y=opt(1)*sinc(2*opt(1).*t)-opt(2)*sinc(2*opt(2).*t);case {'sinc squared'} %generates sinc squared functiony=sinc(2*pi*edge*16.*t/(5*T)).^2;case {'sweep'} %generate frequency sweeptheta=(opt(1)+(opt(2)-opt(1))/T).*(t+T/2);y=real(exp(j*(theta.*(t+T/2)-pi/2)));case {'windowed sweep'} %generate a windows frequency sweeptheta=(opt(1)+(opt(2)-opt(1))/T).*(t+T/2);y=real(exp(j*(theta.*(t+T/2)-pi/2)));c=length(y);edge=min(1,edge);edge=max(0,edge);w=hamming(ceil(c*(1-edge)));w2=[w(1:ceil(length(w)/2));ones(c-length(w),1);w(ceil(length(w)/2)+1:end)]';y=w2.*y;otherwiseerror('invalid pulse type');end%apply a modulationif f~=0y=y.*cos(2*pi*t*f);end%normalize the peak of the pulse to 1p=y./max(abs(y));。

11.微分方程的定义对于duffing 方程032=++x x x ω，先将方程写作⎩⎨⎧--==3112221x x x x x ω function dy=duffing(t,x) omega=1;%定义参数 f1=x(2);f2=-omega^2*x(1)-x(1)^3; dy=[f1;f2];2.微分方程的求解function solve (tstop) tstop=500;%定义时间长度 y0=[0.01;0];%定义初始条件[t,y]=ode45('duffing',tstop,y0,[]);function solve (tstop) step=0.01;%定义步长y0=rand(1,2);%随机初始条件tspan=[0:step:500];%定义时间范围 [t,y]=ode45('duffing',tspan,y0);3.时间历程的绘制时间历程横轴为t ，纵轴为y ，绘制时只取稳态部分。

plot(t,y(:,1));%绘制y 的时间历程 xlabel('t')%横轴为t ylabel('y')%纵轴为y grid;%显示网格线2axis([460 500 -Inf Inf])%图形显示范围设置4.相图的绘制相图的横轴为y ，纵轴为dy/dt ，绘制时也只取稳态部分。

红色部分表示只取最后1000个点。

plot(y(end-1000:end ,1),y(end-1000:end ,2));%绘制y 的时间历程xlabel('y')%横轴为yylabel('dy/dt')%纵轴为dy/dt grid;%显示网格线5.Poincare 映射的绘制对于不同的系统，Poincare 截面的选取方法也不同对于自治系统一般每过其对应线性系统的固有周期，截取一次 对于非自治系统，一般每过其激励的周期，截取一次例程：duffing 方程032=++x x x ω的poincare 映射 function poincare(tstop)global omega; omega=1;T=2*pi/omega;%线性系统的周期或激励的周期 step=T/100;%定义步长为T/100 y0=[0.01;0];%初始条件tspan=[0:step:100*T];%定义时间范围 [t,y]=ode45('duffing',tspan,y0);for i=5000:100:10000%稳态过程每个周期取一个点 plot(y(i,1),y(i,2),'b.'); hold on;% 保留上一次的图形 endxlabel('y');ylabel('dy/dt');3Poincare 映射也可以通过取极值点得到 function poincare(tstop) y0=[0.01;0];tspan=[0:0.01:500];[t,y]=ode45('duffing',tspan,y0); count=find(t>100);%截取稳态过程 y=y(count,:);n=length(y(:,1));%计算点的总数 for i=2:n-1if y(i-1,1)+eps<y(i,1) && y(i,1)>y(i+1,1)+eps % 简单的取出局部最大值plot(y(i,1),y(i,2),'.'); hold on end endxlabel('y');ylabel('dy/dt');6.频谱yy=fft(y(end-1000:end,1)); N=length(yy); power=abs(yy);freq=(1:N-1)*1/step/N;plot(freq(1:N/2),power(1:N/2)); xlabel('f(y)') ylabel('y')7.算例duffing 方程03=++x x x的时间历程，相图，频谱和poincare 映射。

MATLAB弹塑性时程分析法编程弹塑性时程分析是工程结构力学中的一种重要分析方法，用于评估结构在地震等动力荷载下的变形和应力分布。

MATLAB是一种非常强大的科学计算软件，具有丰富的工具箱和函数，可以方便地编写弹塑性时程分析的程序。

本文将介绍如何用MATLAB编程实现弹塑性时程分析。

1.弹塑性分析概述弹塑性分析是一种结构稳定性的计算方法，它考虑了结构的非线性行为，如塑性变形和残余应力。

弹塑性分析的基本思想是将结构划分为弹性和塑性两个部分，根据结构的实际受力情况，逐步计算结构的位移、应力和变形等参数。

2.弹塑性时程分析原理弹塑性时程分析是指以地震动作为输入，计算结构的时程响应。

其基本步骤是：首先，根据结构参数和地震动波特性，建立结构的动力模型。

然后，采用数值积分方法，按照时间步进逐步计算结构的位移、速度和加速度等参数，直到达到要求的计算时间。

在计算过程中，根据结构的非线性本构关系和塑性溃效准则，判断应力状态是否进入塑性阶段，并更新剩余强度等参数。

3.弹塑性时程分析MATLAB编程步骤(1)建立结构的动力模型首先，根据结构的几何形状和材料性质，使用MATLAB建立结构的节点和单元模型。

可以利用网格划分法或几何变换法进行离散化，以获得结构的节点和单元信息。

(2)定义地震动输入根据地震动加速度时程图，使用MATLAB定义地震动输入信号。

可以通过读取实测地震数据，或者使用地震动模拟软件产生地震动波进行模拟。

(3)定义结构的本构关系根据结构的材料性质和截面参数，使用MATLAB定义结构的本构关系。

可以根据结构的线性弹性或非线性塑性材料模型，采用协调变形法或增量处理法进行计算。

(4)制定计算控制策略根据结构的强度要求和计算时间，制定合理的计算控制策略。

这包括选择合适的时间步长和计算时程，以及考虑计算结果的误差控制和稳定性分析。

(5)编写弹塑性时程分析算法根据以上步骤，编写MATLAB程序来实现弹塑性时程分析。

风荷载作用下高层建筑系统的模型分析宿金成，王幼清（哈尔滨工业大学土木工程学院，哈尔滨150090）【摘要】风荷载是随机变化的。

考虑风荷载的随机性以及土－结构相互作用的计算方法对结构进行分析，与实际情况更接近。

文中首先对风荷载的处理方法和土－结构体系的阻尼作分析，然后考虑脉动风的空间相关性应用谐波叠加法编制matlab程序模拟脉动风速时程。

进而，分析刚性地基模型和土－结构模型的动力特性。

以此为基础，对两种模型进行动力和静力分析。

分析表明，地基土对结构整体的动力性能影响显著，土－结构模型的结构顶面最大水平位移比刚性地基模型的大，基底剪应力比刚性地基模型的小。

【关键词】土－结构；脉动风；随机分析；空间相关性【中图分类号】TU973.32【文献标识码】B【文章编号】1001－6864（2012）11－0035－03 HIGH-RISE OVERALL RESPONSE UNDER WIND LOADSU Jin-cheng，WANG You-qing（School of Civil Engineering，Harbin Institute of Technology，Harbin150090，China）Abstract：Wind load changes randomly．Considering the randomness of the wind and soil-structure，structural analysis is realistic．First，wind load and damping of soil-structure were processed．Second，wind speed-time curves were simulated applying harmonic superposition method by matlab．Third，struc-tural dynamic properties were analyzed of rigid foundation and soil-structure model．Finally，static and dynamic analyses of the models were processed．The results show that soil is significantly on structural dynamic property．The maximum horizontal displacement of soil-structure is more than rigid foundation model’s，but base shear stress is less．Key words：soil-structure；wind speed-time；stochastic analysis；spatial correlation风速时程是一个频域较宽且频率变化激烈的时变过程。