PrettyJerry 发表于 2017-9-13 15:11

求解这种图像MATLAB如何画出来的！！！

系统的方程为：


系统的图像为：


图中F的大小分别为：0.002,0.01,0.05,0.1,0.2
求大神解救，真情感谢！！！


linlin820 发表于 2017-9-13 15:45

这个是典型的非线性系统跳跃现象，可以参考下面的文献

PrettyJerry 发表于 2017-9-13 16:50

linlin820 发表于 2017-9-13 15:45
这个是典型的非线性系统跳跃现象，可以参考下面的文献

多谢！

PrettyJerry 发表于 2017-9-17 18:28

linlin820 发表于 2017-9-13 15:45
这个是典型的非线性系统跳跃现象，可以参考下面的文献

你好，文献里面2.1程序我运行出来是这种。和结果差太多。这是什么情况？？？

MVH 发表于 2017-9-18 11:17

PrettyJerry 发表于 2017-9-17 18:28
你好，文献里面2.1程序我运行出来是这种。和结果差太多。这是什么情况？？？

是不是绘图所选取的参数域不全造成的
这种图应该是分段做出来的

PrettyJerry 发表于 2017-9-18 15:06

MVH 发表于 2017-9-18 11:17
是不是绘图所选取的参数域不全造成的
这种图应该是分段做出来的

谢谢帮助，是我编程有问题。修改后已经能画出来了。

陌影 发表于 2017-9-18 16:05

可否把代码分享一下 楼主

PrettyJerry 发表于 2017-9-18 20:56

陌影 发表于 2017-9-18 16:05
可否把代码分享一下 楼主

%强迫振动达芬方程幅频响应曲线的解析求法
%参数设定
epsilon=0.2;%小参数
gammal=4;%非线性参数
F1=0.3;%激励幅值
zetal=0.25;%阻尼项
a=linspace(0.1,2.0,100);%计算响应幅值范围
for ii=1:1:length(a)
    %根据（23）式，对应一个响应幅值，可能存在两个激励频率值；
    %Omega1为较小的激励频率值
    Omega1(ii)=(1+3*epsilon*gammal/8*a(ii)^2-sqrt((epsilon*F1/(2*a(ii)))^2-(epsilon*zetal)^2));
    %Omega1对应的本征函数值
    lamOmega11=sqrt(-((Omega1(ii)-1)/epsilon-3*gammal/8*a(ii)^2)*((Omega1(ii)-1)/epsilon-9*gammal/8*a(ii)^2))-zetal;
    lamOmega12=-sqrt(-((Omega1(ii)-1)/epsilon-3*gammal/8*a(ii)^2)*((Omega1(ii)-1)/epsilon-9*gammal/8*a(ii)^2))-zetal;
    %Omega2为较大的激励频率值
    Omega2(ii)=(1+3*epsilon*gammal/8*a(ii)^2+sqrt((epsilon*F1/(2*a(ii)))^2-(epsilon*zetal)^2));
    %Omega2对应的本征函数值
    lamOmega21=sqrt(-((Omega2(ii)-1)/epsilon-3*gammal/8*a(ii)^2)*((Omega2(ii)-1)/epsilon-9*gammal/8*a(ii)^2))-zetal;
    lamOmega22=-sqrt(-((Omega2(ii)-1)/epsilon-3*gammal/8*a(ii)^2)*((Omega2(ii)-1)/epsilon-9*gammal/8*a(ii)^2))-zetal;
    if lamOmega11==conj(lamOmega21) %如果两个本征值相等退出计算（即此时计算到频响曲线最高点）
   plot(Omega1(ii),a(ii),'k.','MarkerSize',15')
   hold on
   break
    else
      FG=((Omega1(ii)-1)/epsilon-3*gammal/8*a(ii)^2)*((Omega1(ii)-1)/epsilon-9*gammal/8*a(ii)^2)+zetal^2;
      %对应上文中的（31）式
      if FG<0
            plot(Omega1(ii),a(ii),'color',,'MarkerSize',15) %不稳定点
      else
            plot(Omega1(ii),a(ii),'k.','MarkerSize',15) %稳定点
            hold on
      end
      GF=((Omega2(ii)-1)/epsilon-3*gammal/8*a(ii)^2)*((Omega2(ii)-1)/epsilon-9*gammal/8*a(ii)^2)+zetal^2;
      %对应上下文中的（31）式
      if GF<0
            plot(Omega2(ii),a(ii),'.','color',,'MarkerSize',15) %不稳定点
      else
            plot(Omega2(ii),a(ii),'k.','MarkerSize',15)%稳定点
            hold on
      end
      xlabel('\Omega') %x轴标题
      ylabel('a') %y轴标题
      xlim();%x轴范围
      ylim();%y轴范围
    end %满足条件退出循环（计算得到幅频响应曲线顶点）
end %结束对不同a值得计算


sizhiyuan2006 发表于 2017-10-8 16:21

终于解决了，谢谢分享
页: [1]
查看完整版本: 求解这种图像MATLAB如何画出来的！！！