带参数的广义积分数值算法应用的一个实例

solitarytear

移动集中荷载作用下无限长梁稳态振动的挠度解：
  fm=D*beta^4-m*beta^2*v^2+K-i*yita*v*beta;
  W1=p*exp(i*beta*x)/(2*pi*fm);
其中:D,m,K,yita,p为荷载、梁和地基参数，i为虚数，计算时均为已知。
梁的稳态振动的挠度解W(x)为W1关于beta在负无穷大～无穷大上的积分。
    由于被积函数关于beta是振荡衰减的，因此可选择一个固定的积分上下限，将广义积分化为定积分来求得W(x)和x的关系图（此问题东南大学孙璐及浙江大学周华飞等人都有过研究，有结果可参考）。我认为方法有如下三种：1.根据积分定义来求；2.按照本论坛xjzuo版主等提出的方法来求；3.根据Matlab自带的iFFT算法来获得（此法很多文献都提及）。
下面将每种方法的源程序列出：（程序中参数选取和孙璐文一致）
%%%%带参数的广义无穷积分(方法一.按积分的定义求）
     t0=3;
    x=linspace(-t0,t0,121);
    space=0.21;                                
    A=40;
      f=0;
      for beta=-A:space:A              
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 梁、地基、荷载参数
                D=2.3*10^3;                                          
                m=48.2;
                K=6.89*10^7;
                yita=2.3*10^5;
                p=1.05*10^4;
                v=40;               
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                fm=D*beta^4-m*(beta*v)^2+K-i*yita*v*beta;
                fun=1000*p*exp(i*beta*x)*space/(2*pi*fm);   %mm 单位
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                f=f+fun;
        end  
       plot(x,f)
  hold on  

%%%%带参数的广义无穷积分(方法二.xjzuo版主等提出的方法）
            syms beta
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 梁、地基、荷载参数
                D=2.3*10^3;                                          
                m=48.2;
                K=6.89*10^7;
                yita=2.3*10^5;
                p=1.05*10^4;
                v=40;               
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                t0=3;
                x=linspace(-t0,t0,121);
            for j=1:length(x)
                xx=x(j);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                fm=D*beta.^4-m*(beta*v).^2+K-i*yita*v*beta;
                fun=inline(1000*p*exp(i*beta*xx)./(2*pi*fm));   %mm 单位
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
               f(j)=quadl(fun,-40,40);  %% quad实质上即是自适应性simpson积分算法；
           end  
       plot(x,f)
  hold on  

%%%%带参数的广义无穷积分(方法三. 调用Matlab的IFFT算法）
运行后可知，方法一很快得到与孙璐等相同的结果；方法二也可得到满意的结果，但比方法一耗时长，大家可试一试；
至于方法三，可将被积函数积分限截断，并写成级数的形式，将x离散成很多点，此形式和IFFT算法的形式十分相似，许多文献也提及用此法来获得数值结果，编程应该较简单，但我尝试用此思路编程，一直未能得到和方法一，方法二的结果，求高手指点！