请各位对照我师兄的论文看看我的程序

VibrationMaster

还有激振力的频率为25/(约为4Hz)依据是什么?
如果要想"数字"上敏感一些, 应该将它改成系统的固有频率sqrt(K/M)=sqrt(67.51/31.2)=1.471. 这样K减5%肯定会很明显.只不过我不喜欢这样的数学游戏(玩够了),当然你们导师还要让你玩,就玩吧.

猴王

我想我可以称你为老师，但现在的问题是，我不知道下面这个循环语句怎么加进去，麻烦您给我把它的正确位置和其他部分的变化，对照我的程序修正一下谢谢了，您的意思是不是我加到function里面，在执行文件中就不用出现这个循环语句了，但我每次运行的时候，都说k没有附值，但当我把这个语句加进去以后，我就无法识别损伤信号，再就是差分的差别不大，所以我想，您如果有时间，麻烦您能多给我一下指点，无论数字游戏也好，还是科学也好，我现在都不在乎，现在最关键的就是陪导师玩下去
if(t<30)
        K=67.51;
else
        K=67.51*0.95;
end

猴王

当s1=diff(v);s=s1+randn(1,12000),其他的都没有修改，出现了错误，
??? Error using ==> +
Matrix dimensions must agree.

zhlong

s1=diff(v)

由于diff是数值计算 s1比v少了一个数据

VibrationMaster

看看这样行不行？
function xdot=danziyoudukuangjia(t,x)
   
  F=20*sin(1.47*t);
  M=31.2;
  zeta=0.01;

  
  if(t<30)
     K=67.51;
  else
     K=67.51*0.95;
  end
%K=67.51;
  
C=2*zeta*sqrt(K*M);
  
  xdot=[x(2);(1/M)*(F-C*x(2)-K*x(1))];
  
end
==================================
执行文件
%clc;clear all;
close all;
t0=0;
te=60;
dt=0.005;%采样间隔
t=t0:dt:te;

   x0=[0,0]';
   [t,x]=ode45('danziyoudukuangjia',t,x0);
   d=x(:,1);
   v=x(:,2);

s1=(diff(v(1:end-1))+diff(v(2:end)))/2;  %中心差分法
s1=[diff(v(1:2));s1;diff(v(end-1:end))]; %少的两个点用前向和后向差分法补齐
s1=s1/dt; %近似为导数
%T=12001;
%
%x2=zeros(1,T);
%x2(6001)=1;
s=s1+0.01*randn(size(s1));%+x2;

%subplot(311);plot(t,s);
%AXIS([0 60 -10 10]);%单层框架加速度信号输出

%下面进行离散的单尺度小波变换并生成，各尺度上的信号

[c,l]=wavedec(s,5,'bior6.8');%对第一信号进行3尺度一维离散小波分解，采用墨西哥小帽函数
%提取结构的低频和高频信号
ca3=appcoef(c,l,'bior6.8',5);%提取第三尺度系数的低频
[cd1,cd2,cd3,cd4,cd5]=detcoef(c,l,[1,2,3,4,5]);%提取第一、二、三尺度系数的高频
%重构信号的低频和高频部分
a3=wrcoef('a',c,l,'bior6.8',3);
d1=wrcoef('d',c,l,'bior6.8',1);
d2=wrcoef('d',c,l,'bior6.8',2);
d3=wrcoef('d',c,l,'bior6.8',3);
d4=wrcoef('d',c,l,'bior6.8',4);
d5=wrcoef('d',c,l,'bior6.8',5);
%显示多尺度一维信号的分解结果

subplot(611);plot(t,a3);title('信号的低频细节部分');%AXIS([0 60 -5 5]);
subplot(612);plot(t,d1);title('Detailsd1');%AXIS([0 60 -5 5]);
subplot(613);plot(t,d2);title('Detailsd2');%AXIS([0 60 -5 5]);
subplot(614);plot(t,d3);title('Detailsd3');%AXIS([0 60 -5 5]);
subplot(615);plot(t,d4);title('Detailsd4');%AXIS([0 60 -5 5]);
subplot(616);plot(t,d5);title('Detailsd5');%AXIS([0 60 -5 5]);

猴王

十分感谢VibrationMaster老师，谢谢您给我上的这一课，无论怎么样，我真的觉得这段时间跟您和很多人都需到了很多有用的知识，让我受益匪浅，我的师兄当初就是把论文作的不像样子，至于怎么出来的我都不知道，他现在在国外，我联系了，他说很忙也就算了，哎，他是师兄我还能说什么，但导师不管这个什么客观理由，他就要结果，所以我都没办法了，老师，我觉得这个问题的加速度信号还是有些问题，单纯输出加速度的信号时候，开始的振幅，比我师兄的加速度信号低了不少，他的在5左右，而我的也就是一般左右，能不能请您帮帮忙，看看加速度信号该如何处理，是不是要加个初值，也是不是只有这样才能达到这样的效果，还得请你指导阿

VibrationMaster

1. 对你的韧劲表示钦佩. 我现在没有正式的研究生,曾经有几个研究生要问寻求帮助,我的目标都比较难,他们后都不找我了.现在做学术很多都是混混,这不仅在国内,而且国外也是这样,只是比例不一样. 你现在要毕业,看来也得学韦小宝了.如果自己不认真,在所谓的学术圈子里混,其实也挺可怜.
2.对于你的问题,我实在没有什么高招. 如果想自己能重复的话, 而且图形看起来也还算漂亮的话:1)把噪声强度减小;2)尽量用小尺度上的系数;3)小波的滤波器的长度选短一些,比如Harr小波;4)画图的时候,把波动显著尺度的图画的大一些,而不是所有的尺度都用相同大小.
3.按照你的设计,问题的关键是结构损伤后,系统的固有频率发生变化.但是这时响应频率主要取决于系统的激励频率,因此频率应该不会发生显著变化. 速度和位移在ode求解中是要连续的,有可能加速度发生简短.也许加加速度(加速度再求导一次)更明显一些,但是噪声一定要小

猴王

viberationmaster 老师，我现在做的处理是把损伤定义在5秒钟，时间定义在10秒钟，但我的图形出现了比较奇怪的现象，在十秒的时刻有个突变，在第三尺度的细节部分，我想上传我的图片，但格式大于200kb无法上传，所以我把我的改动给您贴出来，请您指教
function xdot=danziyoudukuangjia1(t,x)
   
  F=20*sin(25*t+30);
  M=31.2;
  zeta=0.01;
  
  if(t<5)
     K=67.51;
  else
     K=67.51*0.95;
  end
%K=67.51;
  
C=2*zeta*sqrt(K*M);
  
  xdot=[x(2);(1/M)*(F-C*x(2)-K*x(1))];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%clc;clear all;
close all;
t0=0;
te=10;
dt=0.005;%采样间隔
t=t0:dt:te;
   x0=[0,0]';
   [t,x]=ode45('danziyoudukuangjia',t,x0);
   d=x(:,1);
   v=x(:,2);
s1=(diff(v(1:end-1))+diff(v(2:end)))/2;  %中心差分法
s1=[diff(v(1:2));s1;diff(v(end-1:end))]; %少的两个点用前向和后向差分法补齐
s1=s1/dt; %近似为导数
%T=12001;
%
%x2=zeros(1,T);
%x2(6001)=1;
s=s1+0.01*randn(size(s1));%+x2;
%subplot(311);plot(t,s);
%AXIS([0 60 -10 10]);%单层框架加速度信号输出
%下面进行离散的单尺度小波变换并生成，各尺度上的信号
[c,l]=wavedec(s,5,'bior6.8');%对第一信号进行3尺度一维离散小波分解，采用墨西哥小帽函数
%提取结构的低频和高频信号
ca3=appcoef(c,l,'bior6.8',5);%提取第三尺度系数的低频
[cd1,cd2,cd3,cd4,cd5]=detcoef(c,l,[1,2,3,4,5]);%提取第一、二、三尺度系数的高频
%重构信号的低频和高频部分
a3=wrcoef('a',c,l,'bior6.8',3);
d1=wrcoef('d',c,l,'bior6.8',1);
d2=wrcoef('d',c,l,'bior6.8',2);
d3=wrcoef('d',c,l,'bior6.8',3);
d4=wrcoef('d',c,l,'bior6.8',4);
d5=wrcoef('d',c,l,'bior6.8',5);
%显示多尺度一维信号的分解结果
subplot(611);plot(t,a3);title('信号的低频细节部分');%AXIS([0 60 -5 5]);
subplot(612);plot(t,d1);title('Detailsd1');%AXIS([0 60 -5 5]);
subplot(613);plot(t,d2);title('Detailsd2');%AXIS([0 60 -5 5]);
subplot(614);plot(t,d3);title('Detailsd3');%AXIS([0 60 -5 5]);
subplot(615);plot(t,d4);title('Detailsd4');%AXIS([0 60 -5 5]);
subplot(616);plot(t,d5);title('Detailsd5');%AXIS([0 60 -5 5]);

[ 本帖最后由 zhlong 于 2007-7-8 15:09 编辑 ]

zhlong

花了半天得到的结果，楼主看看下面这些图的结果是否符合你的要求。

[ 本帖最后由 zhlong 于 2007-7-8 14:18 编辑 ]

VibrationMaster

10秒处为边界不连续，与损伤无关。 任何一段数据作变换都会出现这种情况。

猴王

对，就是这样的，希望你能让我看看你的程序，十分感谢

zhlong

你主要是M K C等参数没设置好，可能你没注意它们之间的单位换算。
我修改了一下你的danziyoudukuangjia程序，如下：

F=20*sin(100*t+30);
  M=31.2;
  zeta=0.01;
  
  if(t<5)
     K=67.51*10^5;
  else
     K=67.51*10^5*0.95;
  end
%K=67.51;
  
% C=2*zeta*sqrt(K*M);
C=1.4*10^3;
  
  xdot=[x(2);(1/M)*(F-C*x(2)-K*x(1))];

zhlong

刚才作上面的图后，还没来得及贴程序，死机了。
现在用上面的程序得到的结果有差异，可能还得修改一下。
反正就是你的K和C的数量级太小了，你自己朝着这个方向，不断试试应该就可以得到满意的结果了。

和Vibmaster老师说的一样，这有点玩数字游戏的味道！

[ 本帖最后由 zhlong 于 2007-7-8 19:29 编辑 ]

zhlong

现在用上面的程序得到的结果有差异，可能还得修改一下。

不过你如果

s=s1+0.01*randn(size(s1));%+x2;

这里面不加噪声的话还是能看到突变的。

VibrationMaster

原帖由 zhlong 于 2007-7-8 15:39 发表

不过你如果 这里面不加噪声的话还是能看到突变的。

你所说的突变是不是细节2上的4-5之间的那个东东，它离5太远了。我不知道为什么会出现在这个位置。