Pade approximant in Transport delay block

5 次查看（过去 30 天）

Kashish Pilyal 2022-6-21

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1745045-pade-approximant-in-transport-delay-block

评论： Kashish Pilyal 2022-6-28

I have some queries about the pade approximant in the transport delay block. I know about the pade approximation and how is it done but I am a little confused as to how it is applied here. Also I assume, if in the transport delay block I increase the value of the order of pade approximant, I will get better results at the cost of longer simulation times. Am I assuming right?

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Paul 2022-6-22

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1745045-pade-approximant-in-transport-delay-block#answer_991190

What exactly is the confusion?

Increasing the order of the Pade approximant should have no effect on simulation run time or simulation results. That parameter is only used in developing the linearized model from the simulation block diagram. Have you seen evidence to the contrary?

25 个评论
显示 23更早的评论隐藏 23更早的评论

Paul 2022-6-22

在 MATLAB Online 中打开

Yes, an nth order Pade approximation will add n states to the linearized model

In theory, increasing the order the Pad approximation yields a more accurate approximation of the phase response of the transport delay. But, keep in mind that the Pade approximation is also adding n right-half-plane zeros to your model, which will also have an impact on the results. A 100th order Pade approximation sounds quite excessive and probably will be numerically unreliable, especially if linmod uses a tf representation of the Pade approximation as returned from pade (I don't know that it does and it would be interesting to find out). For example

T = 0.1;

options = bodeoptions;

options.PhaseWrapping = 'on';

options.MagVisible = 'off';

s = tf('s');

h = exp(-s*T);

figure

bodeplot(h,pade(h,2),pade(h,4),logspace(-1,3,500),options)

For a 0.1 second delay, the 2nd order approximant is good to about 30 rad/sec and the 4th order out to about 70-80 rad/sec.

But we run into problems with 100th order approximation

h100 = pade(pade(h,100));
h100.den{:}(end-8:end)
ans = 1×9
1.0e+306 *

    0.0001    0.0118    1.6978       Inf       Inf       Inf       Inf       Inf       Inf

The denominator can't be computed correctly for this value of T!

The bottom line is that the order of the approximation should be as low as possible while still accurately representing the phase of the delay out to some frequency of interest and not causing other numerical problems in the analysis. That frequency of interest will be determined by location of the delay in the system and the properties of the system and/or the properties of the signals being input to the delay.

Paul 2022-6-23

编辑：Paul 2022-6-23

在 MATLAB Online 中打开

Hmm. I've never thought about this question from a theoretical perspective. I'd like to speculate that, in theory, there should be no upper limit on the order of the Pade approximant where instability is induced, just based on a basic Nyquist criterion argument. But that's just speculation on my part.

However, in practice, there is likely a practical upper bound on the order of the approximant before it becomes impractical to use form a computational perspective. As shown above, for T = 0.1, Matlab can't numerically compute the denominator coefficients for N = 100. But problems may arise even for smaller values of N.

T = 0.1;

s = tf('s');

h = exp(-s*T);

pzplot(pade(h,2),pade(h,4),pade(h,20),pade(h,30),ss(pade(h,30)))

legend('N=2','N=4',"N=20","N=30",'N=30 ss')

This plot is attempting to show that the Pade approximant has a pole/zero pattern and that pattern starts to break down somewhere between N = 20 and N = 30 for T = 0.1. Furthermore, the tf form for N=30 can't be reliably converted to ss form, which also indicates concern. As N increases, the poles (and zeros) are migrating further away from origin, the effect of which needs to be understood in the context of any specific analysis.

I'm afraid I can't offer much more on how to determine the upper bound on the order. As stated above, my approach is to use the lowest order such that the phase of the approximant is accurate enough over the frequency range of interest.

Paul 2022-6-26

An H-inf norm greater than 1 might be undesirable for your application, but it does not indicate instability which is why I found several previous comments above very confusing.

This comment stated: "For the pade order of 2 or 3, the delay can be as big as 40 secs for "h=1.3" while maintaining stability of the system. When you increase the pade order to 20, you get slightly lesser values of delay for which system will be stable (the values for which the value of H infinity norm is equal to 1 and does not go larger than 1)."

Based on where this discussion is now, I think it should have stated: For the pade order of 2 or 3, the delay can be as big as 40 secs for "h=1.3" and keep the H infinity norm of the system less than 1. When you increase the pade order to 20, you get slightly lesser values of delay for which system will have H infinity norm less than 1.

Is my modification to the quoted statement accurate?

If so, I'm still not sure what the issue actually is. Apparently the system with h = 1.3 can tolerate nearly 40 seconds ("as bigs as 40 secs" or "sliightly lesser values") of delay before the H infinity norm exceeds 1. If you're concerned that that result is not sensible, then perhaps you need to reassess whether or not the H-infinity norm is the proper figure-of-merit for this problem.

Paul 2022-6-26

在 MATLAB Online 中打开

Here is my attempt to recreate the results in Matlab for h = 1.3. I also used Simulink w/linearize and saw the same answers.

s = tf('s');

ki=1e-05;

kp=0.68;

kd=10;

tau=0.1;

h = 1.3;

tf_a=1/((tau*s)+1);

tf_b=((h*s)+1)/(s^2)*(kp + kd*s + (ki/s));

cons_a=(1-(tau/h));

pid2 = tf([kd kp ki],[1 0 0 0]);

% 20 second delay

Td = 20;

sys20 = series(parallel(tau/h*exp(-s*Td),pid2),feedback(tf_a,parallel(tf_b,-cons_a)));

% H-inf norm for 5th and 20th order Pade

[hinfnorm(minreal(pade(sys20,5))) hinfnorm(minreal(pade(sys20,20)))]

4 states removed. 4 states removed.

ans = 1×2

1.0000 1.0000

bodemag(sys20,minreal(pade(sys20,5)),minreal(pade(sys20,20)));

4 states removed. 4 states removed.

% 40 second delay

Td = 40;

sys40 = series(parallel(tau/h*exp(-s*Td),pid2),feedback(tf_a,parallel(tf_b,-cons_a)));

% H-inf norm for 5th and 20th order Pade

[hinfnorm(minreal(pade(sys40,5))) hinfnorm(minreal(pade(sys40,20)))]

4 states removed. 4 states removed.

ans = 1×2

1.0000 1.0000

bodemag(sys40,minreal(pade(sys40,5)),minreal(pade(sys40,20)));

4 states removed. 4 states removed.

So this code shows H-inf norm of 1 with Td = 20 or 40 and with Pade order of either 5 or 20. Do these results not match yours?

Paul 2022-6-27

在 MATLAB Online 中打开

The posted code does not appear to match the block diagram in the Simulink model.

1. where is cons_a?

2. where is the time delay?

3. how did m = -8.918 make it into the state space model in the A-matrix? As far as I can tell, -8.8918 doesn't multiply anything the block diagram in the Simulink model. If you want to treat that as a gain on the disturbance, then it should show up in a second column in the B-matrix.

4. why is tau not used?

The result from the ss model is:

h=0.5;
m=-8.918;
tau=0.1;
ki=1e-05;
kp=0.15;
kd=10;
A=[0           0           0         0         0;
   0           0           1         0         0;
   0           0           0         1         0;
   m          -ki         -kp       -kd        0;
   0           0           0        1/h       -1/h];
B=[0;
   0;
   0;
   0;
   1];
C=[0      0         0      -1/h    1/h];
D=[0];
sys=ss(A,B,C,D);
tf(minreal(sys))
4 states removed.

ans =
 
    2
  -----
  s + 2
 
Continuous-time transfer function.

Is that what you expect? Is that also what you get from linearizing the Simulink model?

As best I can tell, that result would only be obtained if cons_a = 0, the transport delay is 0, and tau = h = 0.5 (which would make cons_a = 0), and none of these conditions should be true for the problem at hand, which started out with

tau = 0.1;
h = 1.3;
delay = 20; % or something like that

It's very difficult to provide any assistance if the rules of the game keep changing.

Paul 2022-6-27

At this point, I really don't have any more to contribute. Good luck with your project.

Kashish Pilyal 2022-6-28

Ok, thank you for all the assistance. It is really appreciated.

请先登录，再进行评论。

类别

Simulink Simulink Environment Fundamentals Block Libraries Continuous

在 Help Center 和 File Exchange 中查找有关 Continuous 的更多信息

产品

Simulink

版本

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by

Pade approximant in Transport delay block

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

25 个评论
显示 23更早的评论隐藏 23更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

Pade approximant in Transport delay block

0 个评论 显示 -2更早的评论隐藏 -2更早的评论

采纳的回答

25 个评论 显示 23更早的评论隐藏 23更早的评论

更多回答（0 个）

另请参阅

类别

标签

产品

版本

Community Treasure Hunt

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

25 个评论
显示 23更早的评论隐藏 23更早的评论