[自动控制原理] Matlab求解二阶系统阶跃响应(过阻尼)响应时间ts的ts/T1-T1/T2曲线方法(模拟二阶系统的阶跃响应及其性能分析)

学习“二阶系统的阶跃响应的性能分析”之后,每次用T1/T2求解ts时,都'只能查表获得'(这是老师说的哈),为了更方便地求得更加精确地ts/T1值,编写了此函数,也对二阶系统的阶跃响应有了一些了解(我可能自控原理不太好QAQ)

最初我考虑直接使用solve来求解,直接输入函数,但matlab提示错误“Unable to find explicit solution”,错误代码及错误界面如下,此时y解不出来

 

由自控原理的书可以找到ts/T1和T1/T2满足如下关系:

h(t_{s})=1+\frac{1}{\frac{T_{2}}{T_{1}}-1}e^{-\frac{t_{s}}{T_{1}}}+\frac{1}{\frac{T_{1}}{T_{2}}-1}e^{-\frac{T_1}{T_2}\frac{t_{s}}{T_{1}} }=0.95(其中取5%误差)

在“过阻尼”中,T1,T2和阻尼系数\xi有如下关系:
T_1=\frac{1}{\omega_n}\frac{1}{\xi+\surd{(\xi^2-1)}}

T2=\frac{1}{\omega_n}\frac{1}{\xi^2+\surd{(\xi^2-1)}}

故有

\xi_n=\frac{1}{2\omega_n}(\frac{1}{T_1}+\frac{1}{T_2})

注意!一定不要把欠阻尼和过阻尼的公式搞混,本文针对过阻尼情形,而对于欠阻尼情形,有公式t_s=-\frac{ln0.05+\frac{1}{2}ln(1-\xi ^2) }{\xi \omega_n}\approx \frac{3.5}{\xi \omega_n}

下面是作者尝试将“过阻尼”的方法代入“欠阻尼”公式,并用\xi_n表示出\frac{T1}{T2}得到的错误结果(这差不多荒废了作者的一个下午)

 

后来,也尝试过事先替换值的方法,也失败过很多遍

要注意一点(******matlab中的solve函数不是什么都能解出来的!!!***********)

最终,通过简化方程的方法,将正确的函数曲线画了出来

下面附上我修改了好多遍的函数代码(将这个函数保存为FindResponseTime.m,无输入和输出时,绘制响应曲线,而有输入(输入T1/T2)时,绘制带stem的响应曲线,当有输出(如a=FindResponseTime(3.5))时),计算值并返回

function [result]=FindResponseTime(n)
%输入T1/T2  输出ts/T1
    if nargout==0  %无输出值后,绘制曲线
        syms t
        disp("正在解算并绘制响应函数图像,这可能需要2~3分钟时间.........")
        warning off
        
        x=1.05:0.05:8;
        kexi=(1+x)./(2*sqrt(x));
        coe1=1./(1./x-1);
        coe2=1./(x-1);
        y1=[];
        
        for i=1:1:140
            %----------绘图所需的值矩阵----------
            Eq=coe1(i)*exp(-t)+coe2(i)*exp(x(i))^(-t)==-0.05;
            y1=[y1,solve(Eq,t)];
        end
        
        if nargin==0
        %-------计算图像代码:-----------------------
            hold on
            plot(x,y1,'-b','LineWidth',1.3)
            
            for i=1:1:7
                plot(1+i,y1(20*i),'.r',"MarkerSize",15);
                text(1+i,y1(20*i)+0.2,["\zeta=",num2str(kexi(20*i))]);
            end
            xlabel("T1/T2")
            ylabel("t_s/T_1")
            grid;
            title('二阶系统的过阻尼响应曲线')
            axis([1,8,3,4.8]);
       %--------------------------------
        end
        
        if nargin==1
            disp("请注意:绘制线图时,输入的T2/T1要在1~8之间")
            %-------利用回调计算值-----
            y2=FindResponseTime(n);
            %------------------
            hold on
            plot(x,y1,'-b','LineWidth',1.3)
            for i=1:1:7
                plot(1+i,y1(20*i),'.r',"MarkerSize",15);
                text(1+i,y1(20*i)+0.2,["\zeta=",num2str(kexi(20*i))]);
            end
            stem(n,y2)
            t=(n-1)/0.05
%这里我本想用text和num2str将数据写在旁边,但无法写入,错误代码附在注释中
%text(n,y2-0.05,['t_s/T_1=',num2str(y1(t))])
            
            xlabel("T1/T2")
            ylabel("t_s/T_1")
            grid;
            title('二阶系统的过阻尼响应曲线')
            axis([1,8,3,4.8]);
        end
    end
    
    if nargout==1
    %------------------------
            syms t
            x2=n;
            coff1=1/(1/x2-1);
            coff2=1/(x2-1);
            Eq2=coff1*exp(-t)+coff2*exp(x2)^(-t)==-0.05;
            result=solve(Eq2,t);
    %---------计算输入所相对应的值(输入参数错误时自动提示参数不足)---------
    end
end

这里解决问题的方案是:先单独用x计算出两个e^...前面的系数coe1和coe2,再增加t参数,这样先计算出系数之后,solve的运算量就少了很多()

后面的coff是避免重名的

最后说一下使用中的注意事项

对于FindResponseTime(k) 要求k应当在1-8之间(不然会超出图去)

但是a=FindResponseTime(k),你可以随便输入大于1的(如a=FindResponseTime(9)也没问题)(另外这个比其他两个命令运行快)

加上最终成功的运行截图:

 

 

参考的相关网站:

【二阶系统分析】https://mbd.baidu.com/ma/s/MEpamYvN

【3.3二阶系统时间响应】https://mbd.baidu.com/ma/s/9zrTJtkh

参考书目:

《自动控制原理(第七版)》,胡寿松主编

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

程序菜鸟一只

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值