如何解决simulink控制系统仿真中的代数环

您所在的位置:网站首页 职位状态错误怎么解决 如何解决simulink控制系统仿真中的代数环

如何解决simulink控制系统仿真中的代数环

2024-06-01 01:58| 来源: 网络整理| 查看: 265

目录

 1. 什么是代数环

  2. 如何解决代数环

  3. 多个s函数导致的代数环

  4. 源代码

1. 什么是代数环

         在simulink仿真过程中,当输入信号直接取决于输出信号,同时输出信号也直接取决于输入信号时,由于数字计算的时序性,而出现的由于没有输入无法计算输出,没有输出也无法得到输入的“死循环” ,称之为代数环。

        如下图所示,output = func(input+output)。初始时,由于没有output,所以不能计算func函数;但是为了得到output,又必须要计算func。如此往复就形成了代数环。

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_13,color_FFFFFF,t_70,g_se,x_16

2. 如何解决代数环

——连续模型

        对于连续模型,可以在反馈通道添加Memory模块延时

5b032f0514374f0581dceda54ce3474d.png                   watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

——离散模型

        对于离散模型,可以在反馈通道添加delay模块延时

d9129949b3164807afedfe68c2ff86e2.png                watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

 3. 多个s函数导致的代数环

        对于控制系统的仿真来说,很少会出现代数环的问题。下图这种系统就会出现代数环,但是这种系统本身就不存在(至少控制系统不会出现这种玩意)。

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_11,color_FFFFFF,t_70,g_se,x_16

        当simulink中的模块不足以完成建模时,常用的方法就是使用s函数了。如果只是一个s函数建立的模型,例如下图,也不会出现代数环的问题(如果简单粗暴的用s函数建立上图那种模型当我没说)。

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_17,color_FFFFFF,t_70,g_se,x_16

         如果用到多个s函数,如下图,就会提示代数环的警告了。

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

systemP的模型:(源代码放最后啦)

gif.latex?%5Cbegin%7Bbmatrix%7D%20%5Cdot%7Bx%7D1%28t%29%5C%5C%20%5Cdot%7Bx%7D2%28t%29%20%5Cend%7Bbmatrix%7D%20%3D%20%5Cbegin%7Bbmatrix%7D%200%20%26%201%5C%5C%200%20%26%20-1%20%5Cend%7Bbmatrix%7D%5Cbegin%7Bbmatrix%7D%20x1%28t%29%29%5C%5C%20x2%28t%29%29%20%5Cend%7Bbmatrix%7D+%5Cbegin%7Bbmatrix%7D%200%5C%5C1%20%5Cend%7Bbmatrix%7D%20u%28t%29%20-%20%5Cbegin%7Bbmatrix%7D%200%5C%5C-%5Cfrac%7B1%7D%7B2%7D%20%5Cend%7Bbmatrix%7D%20M%28t%29

gif.latex?y%20%3D%20%5Cbegin%7Bbmatrix%7D%20x1%28t%29%5C%5C%20x2%28t%29%5Cend%7Bbmatrix%7D

systemM的模型:(源代码放最后啦)

gif.latex?M1%20%3D%2019.98%20*%20sign%28n0%20-%20n%29%20*%20%280.005%20*%20sin%28a%29%20+%200.02%20*%20cos%28a%29%29       ①

%20%280.0675%5E2%20-%200.055%5E2%29%29%20*%20%28n0%20-%20n%29   ②

 gif.latex?Mw%3D1.7987%20*%20%28tanh%2832.8317%20*%20w%29-tanh%287.4921%20*%20w%29%29%20+%200.9601%20*%20tanh%286.0375%20*%20w%29%20+%200.0378%20*%20w    ③

gif.latex?Mf%20%3D%200.1%20*%20sin%28w%29%20+%200.32       ④

gif.latex?y%20%3D%20M1+Mmud+Mw+Mf

        对于这种情况如果我们加入Memory模块做延时,虽然能解决代数环的问题,但是如果系统本身不是稳定的,延时环节会让系统产生震荡。

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

         将加入Memory和没加Memory的两种输出的波形图对比:(本身搭建的系统就没有意义,所以不分析,将分析另外一个波形图,以此说明差异)

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

        上面的两个模型没有严密的逻辑分析,紧供分析代数环,所以对于系统震荡的问题不好分析,所以在此插入我另外一个仿真中遇到的两个图(也是systemM和systemP构成的仿真,只是系统严谨)

        下图是加入Memory延时,导致机器学习训练收敛性很差

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

        下图是去掉延时,解决代数环后,通过机器学习后的控制波形

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

         很容易发现,第二个的效果很好啦

        延时对不稳定系统的影响还是很大的。对于多个s函数,导致计算时序的影响,出现代数环的问题,我们可以在s函数上做调整。

         在s函数,flag等于0时,也就是初始化时,我们可以看到

sizes.DirFeedthrough

        的设置,输入是否直接影响输出。将其值设为0(修改连续系统的DirFeedthrough,如果修改纯数学公式的会导致仿真错误,本文中systemP为连续系统,systemM为数学公式);仿真时就不会出现代数环的问题,且波形和有代数环时一样。

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

watermark,type_d3F5LXplbmhlaQ,shadow_50,text_Q1NETiBAU3Vkb1JlYm9vdA==,size_20,color_FFFFFF,t_70,g_se,x_16

        我在做仿真的时候,systemP和systemM可以用simulink提供的模块搭建,但是会比较麻烦,所以采用了s函数,后来也做了完全用simulink提供的模块搭建的仿真,控制效果也和s函数出现的效果类似。

4. 源代码

systemP:

%systemP function [sys,x0,str,ts,simStateCompliance] = systemP(t,x,u,flag) switch flag case 0 [sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes; case 1 sys=mdlDerivatives(t,x,u); case 2 sys=mdlUpdate(t,x,u); case 3 sys=mdlOutputs(t,x,u); case 4 sys=mdlGetTimeOfNextVarHit(t,x,u); case 9 sys=mdlTerminate(t,x,u); otherwise DAStudio.error('Simulink:blocks:unhandledFlag', num2str(flag)); end function [sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes sizes = simsizes; sizes.NumContStates = 2; sizes.NumDiscStates = 0; sizes.NumOutputs = 2; sizes.NumInputs = 2; sizes.DirFeedthrough = 1; sizes.NumSampleTimes = 1; sys = simsizes(sizes); x0 = [0;0]; str = []; ts = [0 0]; simStateCompliance = 'UnknownSimState'; function sys=mdlDerivatives(t,x,u) sys = [0 1;0 -1] * x + [0;1]*u(1) - [0;-1/2]*u(2); function sys=mdlUpdate(t,x,u) sys = []; function sys=mdlOutputs(t,x,u) sys = x; function sys=mdlGetTimeOfNextVarHit(t,x,u) sampleTime = 1; sys = t + sampleTime; function sys=mdlTerminate(t,x,u) sys = [];

systemM:

%systemM function [sys,x0,str,ts,simStateCompliance] = systemM(t,x,u,flag) switch flag case 0 [sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes; case 1 sys=mdlDerivatives(t,x,u); case 2 sys=mdlUpdate(t,x,u); case 3 sys=mdlOutputs(t,x,u); case 4 sys=mdlGetTimeOfNextVarHit(t,x,u); case 9 sys=mdlTerminate(t,x,u); otherwise DAStudio.error('Simulink:blocks:unhandledFlag', num2str(flag)); end function [sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes sizes = simsizes; sizes.NumContStates = 0; sizes.NumDiscStates = 0; sizes.NumOutputs = 1; sizes.NumInputs = 3; sizes.DirFeedthrough = 1; sizes.NumSampleTimes = 1; sys = simsizes(sizes); x0 = []; str = []; ts = [0 0]; %[0 0] simStateCompliance = 'UnknownSimState'; function sys=mdlDerivatives(t,x,u) sys = []; function sys=mdlUpdate(t,x,u) sys = []; function sys=mdlOutputs(t,x,u) n0 = u(1); %r/min a = u(2)*pi/180; %degree w = u(3); %rad/s n = w*60/(2*pi); %r/min M1 = 19.98 * sign(n0 - n) * (0.005 * sin(a) + 0.02 * cos(a)); Mmud = 2 / 15 * pi^2 * 0.05 * 1 * (0.0675^2 * 0.055^2 / (0.0675^2 - 0.055^2)) * (n0 - n); Mw = 1.7987 * (tanh(32.8317 * w)-tanh(7.4921 * w)) + 0.9601 * tanh(6.0375 * w) + 0.0378 * w; Mf = 0.1 * sin(w) + 0.32; ML = M1+Mmud+Mw+Mf; sys = ML; function sys=mdlGetTimeOfNextVarHit(t,x,u) sampleTime = 1; sys = t + sampleTime; function sys=mdlTerminate(t,x,u) sys = [];



【本文地址】


今日新闻


推荐新闻


CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3