首页/文章/ 详情

一文搞懂Matlab时频傅里叶转换——修正后重发

6月前浏览7973


昨天发表的《一文搞懂Matlab时频傅里叶转换》一文中存在严重的理论错误,感谢CAE中学生的热心群友提醒,错误修正的地方将在下文详细说明。

众所周知,振动中会用到响应谱,响应谱分析的主要用途是代替瞬态动力学的时程分析,如下图,将一个复杂的时域信号使用拟合的方法,分离成众多的谐波信号,即正弦或余弦信号,由于频谱图省略了各谐波的相位信息,只记录它们的频率与幅值,所以是近似且快速的计算手段。这个内容在教主的动力学课程中说得很详细了。

严格意义上,响应谱分析中的边界条件是一种响应,因为输出的反应谱是位移(速度、加速度或力)与频率之间的关系,这种关系称为频谱。

1 时频转换的必要性 

很多教程说响应谱不需要工程师自己生成,可以直接查询相关标准。在一些标准试验中的确是这样,比如舰船的海浪冲击可以直接查询GJB150相应文件、地震对建筑物的冲击可以查询国标或欧标的建筑规范标准。    

但是大多数时候,我们要计算的频谱并不是标准手册上有的,而是根据我们自身产品在冲击或振动试验中测得的。一般使用传感器测试产品某处的时间历程响应,即时间与位移(或速度、加速度)的时域信号。如果是在固定约束处,也可直接使用振动台的时域信号,而无需传感器。将这个时域信号转换为频域信号,导入有限元软件中对产品进行仿真计算。如果振动台控制器自带频谱输出,那么都不需要时频转换了。

工作中,甲方也不会给你配一个人专门进行信号处理,往往甲方自己也不清楚边界或加载情况,更有甚者试验都让仿真者参与设计,所以不要想着有人给你个频谱,你拿进ansys/abaqus中咔咔一顿操作交差完事,这么轻松为啥甲方自己不搞呢。

在有限元教材中不会讲Matlab的时频转换,在Matlab的教材中又不会讲振动,这个最重要的环节似乎被各位大佬故意遗忘了,我们本文中要探讨的正是这个中间被忽略的部分。

上面划线这部分存在理论错误。先看周老师的说明:

从严格意义上来说响应谱分析中的边界条件不是传统意义上的载荷而是一种响应因为输入的反应谱是位移、速度、加速度或力与频率之间的关系这种关系被称为频谱。例如:在振动试验台上安装4个单自由度的弹簧质量系统,频率分别为f1、f2、f3、f4,且f1<f2<f3<f4。当振动台以频率f1激振时,可得4个系统的位移响应图,再叠加f3频率激振效果,可得4个系统的位移响应图,最后叠加所有频率的综合激振效果,得到的曲线为频谱,如图 4-1-2 所示。  

个人理解,振动台只是起到扫频作用,比如从0~1000HZ扫频,在某些频率下,产品会发生共振而导致振幅(位移、速度或加速度)变大,记录下0~1000HZ下产品所有的振幅,形成响应谱图。  

所以可见,响应谱是产品在连续频率的输入载荷下表现出的自身响应特性。同样0~1000HZ扫频,产品A与产品B的响应谱可能完全不同。

2 Matlab快速傅里叶变换

在Matlab中,时频转换是通过快速傅里叶变换FFT完成了,我们不需要详细研究FFT内部转换过程,只需要直接调用就行了,以下结合实例介绍如何转换。

本文使用的是16版本的Matlab,按照Matlab桀骜不驯的性格,每年升级都会悄悄咪 咪地改动几个命令,所以不能保证适用于其他版本的Matlab。

Step1 数据分析    

如图,打开ZXS.csv文件,可以看到这个信号共1s,每个采样数据间隔0.01s,即采样率100。这里一共有100个数据,注意数据的数量必须是偶数个,数据最上方不能有表头,否则需要修改后文的代码。

Step2 Matlab转换。

首先将工作目录设置到ZXS.csv所在的文件夹。

输入如下代码,%之后的文字为注释。

需要注意的是,需要核对输入文件名称ZXS.csv是否与我们保存的数据文件名相同,核对采样率是否正确。

clc;clear all;close all %清理屏幕与缓存

data = xlsread('ZXS.csv'); %19版本读取命令为data = readmatrix('ZXS.csv')    

t = data(:,1);      %对应第一列

y1 = data(:,2);     %对应第二列

figure

y=y1;         % 读取时域数据

Fs=100;     % 采集频率,1秒采集了多少个点(需要自行根据数据更改)

T=1/Fs;       % 采集时间间隔

N=length(y);  % 采集信号的长度

t=(0:1:N-1)*T;   % 定义整个采集时间点

t=t';            % 转置成列向量

% fft变换

Y=fft(y);          % Y为fft变换结果,复数向量

Y=Y(1:N/2+1);      % 只看变换结果的一半即可

A=abs(Y);          % 复数的幅值(模)

f=(0:1:N/2)*Fs/N;  % 生成频率范围

f=f';              % 转置成列向量

% 幅值修正

A_adj=zeros(N/2+1,1);

A_adj(1)=A(1)/N;      % 频率为0的位置

A_adj(end)=A(end)/N;   % 频率为Fs/2的位置

A_adj(2:end-1)=2*A(2:end-1)/N;

%写出数据    

dataa(:,1)=f;               %频率

dataa(:,2)=A_adj;           %幅值

str = ['DXS.csv'];     %写出文件名称

xlswrite(str,dataa);    %19版本为  writematrix(str,dataa)

% 绘制频率幅值图和频谱图

subplot(2,1,1)

plot(t,y,'LineWidth',1)

xlabel('时间')

ylabel('信号值')

title('时域信号')

subplot(2,1,2)

plot(f,A_adj,'LineWidth',3)

xlabel('频率(Hz)')

ylabel('幅值(修正后)')

title('FFT变换幅值图')

输入后,点击上方运行工具,将自动导出转换的频域信号DXS.csv,并生成图像,上面显示的是输入的时域图,下面显示的是输出的频域图。可以看到,这个型号是由三个信号组成,分别是频率=1幅值=5的信号1,频率=8幅值=1的信号2,频率=10幅值=0.5的信号3。   

到输出的DXS.csv文件中,插入散点图能更清晰看到组成。

可以看到,频域信号底部前后各扩展了1hz,形成了三角形,这是因为频谱分辨率导致,分辨率=1/总时间,我们的信号总时间是1s,所以分辨率是1Hz。

现在,你可以将自己的时频信号导入测试

动力学的仿真实操其实很简单,与静力学无异,难就难在理论的解读,我也感觉到自身理论的短板,也请大家多大指教批评。

本文特别感谢卧虎藏龙的CAE中学生的初一(2)班的jiff、黄荣等大佬的热心解答。

参考书籍:《ANSYS  Workbench  有限元分析实例详解  动力学》—周炬、苏金英

参考视频:《ANSYS高级视频教程》——张晔


来源:CAE中学生
WorkbenchAbaqus静力学瞬态动力学振动建筑MATLAB理论控制试验ANSYS
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2024-05-10
最近编辑:6月前
CAE无剑
硕士 | 仿真工程师 CAE中学生
获赞 688粉丝 1507文章 250课程 0
点赞
收藏
作者推荐

10 2D如何划分结构化网格总结——HyperMesh速通10

本文不具备实战意义,主要在于探讨结构化网格的一些技术问题。HM程序越来越聪明,并不需要用户投入巨大的精力,特别适合于LSDYNA、abaqus等软件的前处理。虽然HM程序已经很聪明了,但是很多时候还是需要工程师掌握一定的技巧,才能划分出结构化网格。下面讨论几个网格划分的问题。1 边与点对网格的影响 当面中存在共享边、自由边、T形边等时,对网格划分会有一定的影响。打开HM软件,选择默认的HM环境。打开edge.hm文件。通过by topo观察。使用2D——automesh,选择所有面,使用10mm网格划分尺寸,所有选项默认。可以看到,压缩边不影响网格划分,共享边和T形边上会生成网格节点,自由边也能生成节点,而且两边可以调节不同的份数,自由边上的节点没有共享拓扑。新建Model,打开point.hm文件。打开硬点显示工具,可以看到左侧面含有若干临时节点,显示为黄色大圆点;中间面含有若干自由点,向思维x;右边面含有若干硬点,显示为小圆点。 使用2D——automesh,选择所有面,使用10mm网格划分尺寸,所有选项默认。 可以看到临时节点与自由点都不影响网格的划分,而硬点会强制关联到网格的节点。2 不同形状的映射网格 标准的圆形、等边三角形、四边形都能通过automesh自动映射为结构网格。椭圆、不规则三角形、不规则四边形、五边形可以通过automesh的二级面板调整为全结构网格。注意,各边网格份数需要设置为偶数,才能更容易划分全四边形网格。 更不规则的形状需要进行切分,切分的原则是每一部分有4个边,比如圆的O形切分,三角形的Y形切分等。内圆外方一般采用1/8切割。如果模型中只有半个圆或1/4圆,则程序无法自动映射,也需要切割。对于面上的圆孔,如果不是特别要求,也可以不切分,而是采用quick edit——washer扩展一圈,形成1~3层孔边网格,这一思路是借鉴了流体网格的划分方法。3 变分网格 复杂模型中不可避免地需要在粗网格与细网格之间过渡,如何来保证过渡网格结构规则呢?这就需要用到变分网格,变分网格推荐1分2,1分3等结构。由于1分2的过渡单元只有一个角为锐角,而1分3的过渡单元有两个角为锐角,所以更推荐1分2的变分,还可以进行多次变分。 4 实例详解 最后以一个实例来进行结构网格划分,这个实例是在官方模型,但是每个工程师的划分方法都有点差别,所以本例仅作参考。Step1 选择起始面使用HM默认环境打开midsurface1.hm。我们从4个孔这一方开始划分,使用quick edit——washer将孔扩展一圈1mm区域。 为了更便于控制,将面如下分割。使用automesh对附近面直接划分1mm网格。在二级面板中调节网格份数,对圆边划分12等分,其余边上的节点数量一一对应。除了红圈内,其余网格都较规整。要使红圈内网格规整有很多办法,我们本文使用二级面板中的偏置方法,稍微调整这两条边的偏置系数,点击鼠标中键确定重新划分。 对于强迫症来说有些节点排列不够整齐,我们可以通过Geom——node edit——align node将它们排列整齐。Step2 变分与镜像接下来演示下变分网格,使用两侧斜面来做过渡面,将面按每两个网格切分一次,使用的工具是Geom——quick edit——split surf-line。对这些分割面划分网格。对边都网格数量均设置为1,侧边依此按1-2-1-2这样的规律设置。 对面斜边也可以这样设置,但是由于两个斜边是完全对称的,我们直接镜像过去。工具选择Tool——reflect,选择斜边上的变分网格,再次点击elems框选择duplicate复制到本组件,选择方向为3点定义,N1与N2点选择几何边定义方向,对称基点B选择圆孔切分点。此时镜像过来的网格节点并未与几何绑定,需要通过Geom——node edit——associate绑定节点与几何。是不是这样就完成了?那又不是。节点虽然与几何绑定了,但是与原有的网格并未合并,什么意思呢?放大平面与斜面网格处,可以看到节点并未重合。 需要通过tool——edges合并节点,这一步比绑定几何更重要。注意容差设置不能超过最小单元尺寸,但是要大于节点错开的距离。合并节点后共享边上的节点便全部被一一合并。是不是复制移动或镜像后,共享边的节点全部一一重合就不需要合并了呢?当然不是,这时候节点虽然重合,但是没有合并,在同一位置有两个节点。Step3 划分大面粗网格划分边沿面粗网格。一些转角处的细碎边使用合并节点消除,大面一边切割一边划分网格。Step4 曲面网格细化在悬臂面曲率处细化网格,还是使用变分网格的老套路。将圆孔附近切割出来,切割的原则是使每一边尽量相等 然后将网格从曲率处依次划分过来最后划分圆孔附近的网格。至此网格划分完成。最终网格如下图。模型已经放入这一文的附件中,网格检查、厚度属性赋予等留到下一文实例介绍。 注:使用quick edit、surfs edit等工具切割曲面时,经常会遇到切割跑偏的问题,此时不用纠结,使用Ctrl+Z返回上一步,将模型局部放大,换个姿势,重复一次操作可能就正常了。如果还是跑偏,使用替换节点等方法来处理。网格划分过程中还要注意随时保存,这都是教训。 来源:CAE中学生

未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈